Function to load data

load_data<-function(){
  df<-read.csv("data2/data_clean.csv")
  df$Time.upper.G1 <- ifelse(df$Time.lower.G1 == df$Time.upper.G1, Inf, df$Time.upper.G1)
  df$Time.upper.G2 <- ifelse(df$Time.lower.G2 == df$Time.upper.G2, Inf, df$Time.upper.G2)
  df$Time.upper.G3 <- ifelse(df$Time.lower.G3 == df$Time.upper.G3, Inf, df$Time.upper.G3)
  df$Transplant.Date<-anydate(df$Transplant.Date)
  df$Date.of.Graft.failure<-anydate(df$Date.of.Graft.failure)
  df$DOB<-anydate(df$DOB)
  df$age<-as.numeric(difftime(df$Transplant.Date, df$DOB, units="days"))/365.25
  df<-df %>% filter(!df$Pre.DSAs=="Yes")
  df$time2fail<-ifelse(is.na(df$Date.of.Graft.failure), anydate("2024-03-10")-df$Transplant.Date, df$Date.of.Graft.failure-df$Transplant.Date)
  df$status<-ifelse(is.na(df$Date.of.Graft.failure), 0, 1)
  df$tot_mismatch<-df$ClassI....Mismatch..All+df$ClassII....Mismatch..All
  df$Donor.TGLN<-df$Donor.TGLN %>% as.character()
  df$MRN<-df$MRN %>% as.character()
  df$trans_id<-paste0(df$Donor.TGLN, df$MRN) %>% as.numeric()
  df <- df %>% mutate(row_number = row_number())
  df$DRB<-df$DRB1.Mismatch.No.+df$DRB345.Mismatch.No.
  return(df)
}

Independant variables in question for analysis are saved. So are HLA mismatch variables

HLA_vars<-c("A.Mismatch.No.", "B.Mismatch.No.", "C.Mismatch.No.", "DRB1.Mismatch.No.", "DRB345.Mismatch.No.", "DP.Mismatch.No.", "DQ.Mismatch.No.")
independent_vars<-c(HLA_vars, "age")
cumulative_cols<-c("ClassII....Mismatch..All", "ClassI....Mismatch..All", "tot_mismatch")

HLA_al_vars<-c("A", "B", "C", "DRB1", "DRB345", "DP", "DQ")
independent_al_vars<-c(HLA_al_vars, "age")

load data and plot Kaplan-Meier curve

df<-load_data()
## Warning in paste0(df$Donor.TGLN, df$MRN) %>% as.numeric(): NAs introduced by
## coercion
plot(survfit(Surv(time2fail, status) ~ 1, data = df), xlab = "Time (days)", ylab = "Survival Probability", main = "Kaplan-Meier Estimate")

Subsetting df to make a correlation plot to check for possible confounding.

indep_df<-df[independent_vars]
indep_df%>%cor()%>%corrplot()

svg("allele_cor.svg")
indep_df_al<-df[independent_al_vars]
indep_df_al%>%cor()%>%corrplot()
dev.off()
## png 
##   2

Log-rank tests and Kaplan-Meier curves

plot_quartiles_and_logrank <- function(data, independent_vars, status_var, time_var) {
  results <- list()
  
  for (var in independent_vars) {
    data[var]<-ntile(data[var], 2)
    fit_formula <- as.formula(paste0("Surv(", time_var, ", ", status_var, ") ~ ", var))
    fit <- surv_fit(fit_formula, data = data)
    p <- ggsurvplot(
      fit, 
      data = data, 
      xlab = "Time", 
      ylab = "Survival probability", 
      palette = c("#4DAD5B","#FC4E07"),
      conf.int = TRUE,
      surv.scale = "percent",
    )
    
    log_rank_result <- survdiff(fit_formula, data = data)

    # Store results
    results[[var]] <- list(plot = p, log_rank_result = log_rank_result)
    cat("Log-rank test result for", var, ":\n")
    print(log_rank_result)
    print(p)
    #ggsave(paste0("figs/km_curve/", var, ".svg"))
  }

  return(results)
}

Kaplan Meier curves and log-rank tests are performed for eplet mismatches

result<-plot_quartiles_and_logrank(df, independent_vars, "status", "time2fail")
## Log-rank test result for A.Mismatch.No. :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
##                   N Observed Expected (O-E)^2/E (O-E)^2/V
## A.Mismatch.No.=1 49        2     4.72      1.57      3.31
## A.Mismatch.No.=2 49        7     4.28      1.72      3.31
## 
##  Chisq= 3.3  on 1 degrees of freedom, p= 0.07

## Log-rank test result for B.Mismatch.No. :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
##                   N Observed Expected (O-E)^2/E (O-E)^2/V
## B.Mismatch.No.=1 49        3     4.37     0.432     0.841
## B.Mismatch.No.=2 49        6     4.63     0.409     0.841
## 
##  Chisq= 0.8  on 1 degrees of freedom, p= 0.4

## Log-rank test result for C.Mismatch.No. :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
##                   N Observed Expected (O-E)^2/E (O-E)^2/V
## C.Mismatch.No.=1 49        3     4.39     0.439     0.857
## C.Mismatch.No.=2 49        6     4.61     0.418     0.857
## 
##  Chisq= 0.9  on 1 degrees of freedom, p= 0.4

## Log-rank test result for DRB1.Mismatch.No. :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
##                      N Observed Expected (O-E)^2/E (O-E)^2/V
## DRB1.Mismatch.No.=1 49        4      4.5    0.0549      0.11
## DRB1.Mismatch.No.=2 49        5      4.5    0.0549      0.11
## 
##  Chisq= 0.1  on 1 degrees of freedom, p= 0.7

## Log-rank test result for DRB345.Mismatch.No. :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
##                        N Observed Expected (O-E)^2/E (O-E)^2/V
## DRB345.Mismatch.No.=1 49        6     4.44     0.550      1.09
## DRB345.Mismatch.No.=2 49        3     4.56     0.535      1.09
## 
##  Chisq= 1.1  on 1 degrees of freedom, p= 0.3

## Log-rank test result for DP.Mismatch.No. :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
##                    N Observed Expected (O-E)^2/E (O-E)^2/V
## DP.Mismatch.No.=1 49        4     4.38    0.0337    0.0659
## DP.Mismatch.No.=2 49        5     4.62    0.0321    0.0659
## 
##  Chisq= 0.1  on 1 degrees of freedom, p= 0.8

## Log-rank test result for DQ.Mismatch.No. :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
##                    N Observed Expected (O-E)^2/E (O-E)^2/V
## DQ.Mismatch.No.=1 49        6     4.55     0.464      0.94
## DQ.Mismatch.No.=2 49        3     4.45     0.474      0.94
## 
##  Chisq= 0.9  on 1 degrees of freedom, p= 0.3

## Log-rank test result for age :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
##        N Observed Expected (O-E)^2/E (O-E)^2/V
## age=1 49        4     4.49    0.0533     0.106
## age=2 49        5     4.51    0.0530     0.106
## 
##  Chisq= 0.1  on 1 degrees of freedom, p= 0.7

result2<-plot_quartiles_and_logrank(df, cumulative_cols, "status", "time2fail")
## Log-rank test result for ClassII....Mismatch..All :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
##                             N Observed Expected (O-E)^2/E (O-E)^2/V
## ClassII....Mismatch..All=1 49        7     4.43      1.49      2.94
## ClassII....Mismatch..All=2 49        2     4.57      1.45      2.94
## 
##  Chisq= 2.9  on 1 degrees of freedom, p= 0.09

## Log-rank test result for ClassI....Mismatch..All :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
##                            N Observed Expected (O-E)^2/E (O-E)^2/V
## ClassI....Mismatch..All=1 49        2     4.49      1.38      2.76
## ClassI....Mismatch..All=2 49        7     4.51      1.38      2.76
## 
##  Chisq= 2.8  on 1 degrees of freedom, p= 0.1

## Log-rank test result for tot_mismatch :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
##                 N Observed Expected (O-E)^2/E (O-E)^2/V
## tot_mismatch=1 49        5     4.39    0.0856     0.167
## tot_mismatch=2 49        4     4.61    0.0814     0.167
## 
##  Chisq= 0.2  on 1 degrees of freedom, p= 0.7

Kaplan Meier curves and log-rank tests are performed for allele mismatches

result<-plot_quartiles_and_logrank(df, independent_al_vars, "status", "time2fail")
## Log-rank test result for A :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
##      N Observed Expected (O-E)^2/E (O-E)^2/V
## A=1 49        6     4.18     0.796      1.49
## A=2 49        3     4.82     0.689      1.49
## 
##  Chisq= 1.5  on 1 degrees of freedom, p= 0.2

## Log-rank test result for B :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
##      N Observed Expected (O-E)^2/E (O-E)^2/V
## B=1 49        6     3.94      1.08      1.98
## B=2 49        3     5.06      0.84      1.98
## 
##  Chisq= 2  on 1 degrees of freedom, p= 0.2

## Log-rank test result for C :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
##      N Observed Expected (O-E)^2/E (O-E)^2/V
## C=1 49        3     4.37     0.432      0.84
## C=2 49        6     4.63     0.408      0.84
## 
##  Chisq= 0.8  on 1 degrees of freedom, p= 0.4

## Log-rank test result for DRB1 :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
##         N Observed Expected (O-E)^2/E (O-E)^2/V
## DRB1=1 49        4     4.35    0.0280    0.0543
## DRB1=2 49        5     4.65    0.0262    0.0543
## 
##  Chisq= 0.1  on 1 degrees of freedom, p= 0.8

## Log-rank test result for DRB345 :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
##           N Observed Expected (O-E)^2/E (O-E)^2/V
## DRB345=1 49        4     4.73     0.112     0.236
## DRB345=2 49        5     4.27     0.124     0.236
## 
##  Chisq= 0.2  on 1 degrees of freedom, p= 0.6

## Log-rank test result for DP :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
##       N Observed Expected (O-E)^2/E (O-E)^2/V
## DP=1 49        4     4.12   0.00370    0.0069
## DP=2 49        5     4.88   0.00313    0.0069
## 
##  Chisq= 0  on 1 degrees of freedom, p= 0.9

## Log-rank test result for DQ :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
##       N Observed Expected (O-E)^2/E (O-E)^2/V
## DQ=1 49        6     4.55     0.459     0.932
## DQ=2 49        3     4.45     0.470     0.932
## 
##  Chisq= 0.9  on 1 degrees of freedom, p= 0.3

## Log-rank test result for age :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
##        N Observed Expected (O-E)^2/E (O-E)^2/V
## age=1 49        4     4.49    0.0533     0.106
## age=2 49        5     4.51    0.0530     0.106
## 
##  Chisq= 0.1  on 1 degrees of freedom, p= 0.7

Cox proportional hazards univariate modeling function

univariate_coxph <- function(data, independent_vars, status_var, time_var) {
  results <- list()
  summary_table <- data.frame(Variable = character(),
                              HR = numeric(),
                              Lower_CI = numeric(),
                              Upper_CI = numeric(),
                              p_value = numeric(),
                              C_index = numeric(),
                              stringsAsFactors = FALSE)
  
  for (var in independent_vars) {
    fit_formula <- as.formula(paste0("Surv(", time_var, ", ", status_var, ") ~ ", var))
    fit <- coxph(fit_formula, data = data, na.action = na.exclude)  # Exclude NA values
    p<-ggcoxzph(cox.zph(fit))
    print(p)
    # Extract summary information
    cox_summary <- summary(fit)
    hr <- exp(coef(fit))
    ci <- exp(confint(fit))
    p_val <- cox_summary$coefficients[,5]
    c_index <- cox_summary$concordance[1]  # Extract C-index
    
    # Print model summary for each variable
    cat("\nCox proportional hazards model for", var, ":\n")
    print(cox_summary)
    
    # Store results in summary table
    summary_table <- rbind(summary_table, data.frame(Variable = var, 
                                                     HR = hr, 
                                                     Lower_CI = ci[1], 
                                                     Upper_CI = ci[2], 
                                                     p_value = p_val,
                                                     C_index = c_index))
    
    results[[var]] <- list(Model = fit, Summary = cox_summary)
  }
  
  # Print summary table of key results
  print(summary_table)
  
  return(results)
}

Univariate coxph models for eplet time to failure

cox_results<-univariate_coxph(df, independent_vars, "status", "time2fail")

## 
## Cox proportional hazards model for A.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##                   coef exp(coef) se(coef)    z Pr(>|z|)
## A.Mismatch.No. 0.06325   1.06529  0.04828 1.31     0.19
## 
##                exp(coef) exp(-coef) lower .95 upper .95
## A.Mismatch.No.     1.065     0.9387    0.9691     1.171
## 
## Concordance= 0.682  (se = 0.1 )
## Likelihood ratio test= 1.69  on 1 df,   p=0.2
## Wald test            = 1.72  on 1 df,   p=0.2
## Score (logrank) test = 1.71  on 1 df,   p=0.2

## 
## Cox proportional hazards model for B.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##                   coef exp(coef) se(coef)    z Pr(>|z|)
## B.Mismatch.No. 0.01369   1.01378  0.05469 0.25    0.802
## 
##                exp(coef) exp(-coef) lower .95 upper .95
## B.Mismatch.No.     1.014     0.9864    0.9107     1.128
## 
## Concordance= 0.525  (se = 0.096 )
## Likelihood ratio test= 0.06  on 1 df,   p=0.8
## Wald test            = 0.06  on 1 df,   p=0.8
## Score (logrank) test = 0.06  on 1 df,   p=0.8

## 
## Cox proportional hazards model for C.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##                   coef exp(coef) se(coef)     z Pr(>|z|)
## C.Mismatch.No. 0.09009   1.09427  0.05779 1.559    0.119
## 
##                exp(coef) exp(-coef) lower .95 upper .95
## C.Mismatch.No.     1.094     0.9139    0.9771     1.226
## 
## Concordance= 0.656  (se = 0.076 )
## Likelihood ratio test= 2.28  on 1 df,   p=0.1
## Wald test            = 2.43  on 1 df,   p=0.1
## Score (logrank) test = 2.44  on 1 df,   p=0.1

## 
## Cox proportional hazards model for DRB1.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##                        coef exp(coef)  se(coef)     z Pr(>|z|)
## DRB1.Mismatch.No. 0.0006791 1.0006793 0.0527482 0.013     0.99
## 
##                   exp(coef) exp(-coef) lower .95 upper .95
## DRB1.Mismatch.No.     1.001     0.9993    0.9024      1.11
## 
## Concordance= 0.549  (se = 0.078 )
## Likelihood ratio test= 0  on 1 df,   p=1
## Wald test            = 0  on 1 df,   p=1
## Score (logrank) test = 0  on 1 df,   p=1

## 
## Cox proportional hazards model for DRB345.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##                         coef exp(coef) se(coef)      z Pr(>|z|)
## DRB345.Mismatch.No. -0.06852   0.93378  0.04904 -1.397    0.162
## 
##                     exp(coef) exp(-coef) lower .95 upper .95
## DRB345.Mismatch.No.    0.9338      1.071    0.8482     1.028
## 
## Concordance= 0.593  (se = 0.09 )
## Likelihood ratio test= 2.25  on 1 df,   p=0.1
## Wald test            = 1.95  on 1 df,   p=0.2
## Score (logrank) test = 1.98  on 1 df,   p=0.2

## 
## Cox proportional hazards model for DP.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##                      coef exp(coef)  se(coef)      z Pr(>|z|)
## DP.Mismatch.No. -0.001622  0.998379  0.046560 -0.035    0.972
## 
##                 exp(coef) exp(-coef) lower .95 upper .95
## DP.Mismatch.No.    0.9984      1.002    0.9113     1.094
## 
## Concordance= 0.544  (se = 0.096 )
## Likelihood ratio test= 0  on 1 df,   p=1
## Wald test            = 0  on 1 df,   p=1
## Score (logrank) test = 0  on 1 df,   p=1

## 
## Cox proportional hazards model for DQ.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##                     coef exp(coef) se(coef)      z Pr(>|z|)
## DQ.Mismatch.No. -0.02535   0.97497  0.05212 -0.486    0.627
## 
##                 exp(coef) exp(-coef) lower .95 upper .95
## DQ.Mismatch.No.     0.975      1.026    0.8803      1.08
## 
## Concordance= 0.53  (se = 0.063 )
## Likelihood ratio test= 0.24  on 1 df,   p=0.6
## Wald test            = 0.24  on 1 df,   p=0.6
## Score (logrank) test = 0.24  on 1 df,   p=0.6

## 
## Cox proportional hazards model for age :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##          coef exp(coef)  se(coef)      z Pr(>|z|)
## age -0.002506  0.997497  0.029763 -0.084    0.933
## 
##     exp(coef) exp(-coef) lower .95 upper .95
## age    0.9975      1.003     0.941     1.057
## 
## Concordance= 0.477  (se = 0.111 )
## Likelihood ratio test= 0.01  on 1 df,   p=0.9
## Wald test            = 0.01  on 1 df,   p=0.9
## Score (logrank) test = 0.01  on 1 df,   p=0.9
## 
##                                Variable        HR  Lower_CI Upper_CI   p_value
## A.Mismatch.No.           A.Mismatch.No. 1.0652930 0.9691025 1.171031 0.1902104
## B.Mismatch.No.           B.Mismatch.No. 1.0137822 0.9107412 1.128481 0.8023563
## C.Mismatch.No.           C.Mismatch.No. 1.0942703 0.9770878 1.225507 0.1190252
## DRB1.Mismatch.No.     DRB1.Mismatch.No. 1.0006793 0.9023926 1.109671 0.9897287
## DRB345.Mismatch.No. DRB345.Mismatch.No. 0.9337766 0.8482019 1.027985 0.1623660
## DP.Mismatch.No.         DP.Mismatch.No. 0.9983792 0.9113051 1.093773 0.9722078
## DQ.Mismatch.No.         DQ.Mismatch.No. 0.9749681 0.8802975 1.079820 0.6266642
## age                                 age 0.9974973 0.9409729 1.057417 0.9329038
##                       C_index
## A.Mismatch.No.      0.6818774
## B.Mismatch.No.      0.5247718
## C.Mismatch.No.      0.6558018
## DRB1.Mismatch.No.   0.5488918
## DRB345.Mismatch.No. 0.5932203
## DP.Mismatch.No.     0.5436767
## DQ.Mismatch.No.     0.5299870
## age                 0.4771838
cox_results
## $A.Mismatch.No.
## $A.Mismatch.No.$Model
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##                   coef exp(coef) se(coef)    z    p
## A.Mismatch.No. 0.06325   1.06529  0.04828 1.31 0.19
## 
## Likelihood ratio test=1.69  on 1 df, p=0.1939
## n= 98, number of events= 9 
## 
## $A.Mismatch.No.$Summary
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##                   coef exp(coef) se(coef)    z Pr(>|z|)
## A.Mismatch.No. 0.06325   1.06529  0.04828 1.31     0.19
## 
##                exp(coef) exp(-coef) lower .95 upper .95
## A.Mismatch.No.     1.065     0.9387    0.9691     1.171
## 
## Concordance= 0.682  (se = 0.1 )
## Likelihood ratio test= 1.69  on 1 df,   p=0.2
## Wald test            = 1.72  on 1 df,   p=0.2
## Score (logrank) test = 1.71  on 1 df,   p=0.2
## 
## 
## 
## $B.Mismatch.No.
## $B.Mismatch.No.$Model
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##                   coef exp(coef) se(coef)    z     p
## B.Mismatch.No. 0.01369   1.01378  0.05469 0.25 0.802
## 
## Likelihood ratio test=0.06  on 1 df, p=0.8028
## n= 98, number of events= 9 
## 
## $B.Mismatch.No.$Summary
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##                   coef exp(coef) se(coef)    z Pr(>|z|)
## B.Mismatch.No. 0.01369   1.01378  0.05469 0.25    0.802
## 
##                exp(coef) exp(-coef) lower .95 upper .95
## B.Mismatch.No.     1.014     0.9864    0.9107     1.128
## 
## Concordance= 0.525  (se = 0.096 )
## Likelihood ratio test= 0.06  on 1 df,   p=0.8
## Wald test            = 0.06  on 1 df,   p=0.8
## Score (logrank) test = 0.06  on 1 df,   p=0.8
## 
## 
## 
## $C.Mismatch.No.
## $C.Mismatch.No.$Model
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##                   coef exp(coef) se(coef)     z     p
## C.Mismatch.No. 0.09009   1.09427  0.05779 1.559 0.119
## 
## Likelihood ratio test=2.28  on 1 df, p=0.1308
## n= 98, number of events= 9 
## 
## $C.Mismatch.No.$Summary
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##                   coef exp(coef) se(coef)     z Pr(>|z|)
## C.Mismatch.No. 0.09009   1.09427  0.05779 1.559    0.119
## 
##                exp(coef) exp(-coef) lower .95 upper .95
## C.Mismatch.No.     1.094     0.9139    0.9771     1.226
## 
## Concordance= 0.656  (se = 0.076 )
## Likelihood ratio test= 2.28  on 1 df,   p=0.1
## Wald test            = 2.43  on 1 df,   p=0.1
## Score (logrank) test = 2.44  on 1 df,   p=0.1
## 
## 
## 
## $DRB1.Mismatch.No.
## $DRB1.Mismatch.No.$Model
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##                        coef exp(coef)  se(coef)     z    p
## DRB1.Mismatch.No. 0.0006791 1.0006793 0.0527482 0.013 0.99
## 
## Likelihood ratio test=0  on 1 df, p=0.9897
## n= 98, number of events= 9 
## 
## $DRB1.Mismatch.No.$Summary
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##                        coef exp(coef)  se(coef)     z Pr(>|z|)
## DRB1.Mismatch.No. 0.0006791 1.0006793 0.0527482 0.013     0.99
## 
##                   exp(coef) exp(-coef) lower .95 upper .95
## DRB1.Mismatch.No.     1.001     0.9993    0.9024      1.11
## 
## Concordance= 0.549  (se = 0.078 )
## Likelihood ratio test= 0  on 1 df,   p=1
## Wald test            = 0  on 1 df,   p=1
## Score (logrank) test = 0  on 1 df,   p=1
## 
## 
## 
## $DRB345.Mismatch.No.
## $DRB345.Mismatch.No.$Model
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##                         coef exp(coef) se(coef)      z     p
## DRB345.Mismatch.No. -0.06852   0.93378  0.04904 -1.397 0.162
## 
## Likelihood ratio test=2.25  on 1 df, p=0.1332
## n= 98, number of events= 9 
## 
## $DRB345.Mismatch.No.$Summary
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##                         coef exp(coef) se(coef)      z Pr(>|z|)
## DRB345.Mismatch.No. -0.06852   0.93378  0.04904 -1.397    0.162
## 
##                     exp(coef) exp(-coef) lower .95 upper .95
## DRB345.Mismatch.No.    0.9338      1.071    0.8482     1.028
## 
## Concordance= 0.593  (se = 0.09 )
## Likelihood ratio test= 2.25  on 1 df,   p=0.1
## Wald test            = 1.95  on 1 df,   p=0.2
## Score (logrank) test = 1.98  on 1 df,   p=0.2
## 
## 
## 
## $DP.Mismatch.No.
## $DP.Mismatch.No.$Model
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##                      coef exp(coef)  se(coef)      z     p
## DP.Mismatch.No. -0.001622  0.998379  0.046560 -0.035 0.972
## 
## Likelihood ratio test=0  on 1 df, p=0.9722
## n= 98, number of events= 9 
## 
## $DP.Mismatch.No.$Summary
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##                      coef exp(coef)  se(coef)      z Pr(>|z|)
## DP.Mismatch.No. -0.001622  0.998379  0.046560 -0.035    0.972
## 
##                 exp(coef) exp(-coef) lower .95 upper .95
## DP.Mismatch.No.    0.9984      1.002    0.9113     1.094
## 
## Concordance= 0.544  (se = 0.096 )
## Likelihood ratio test= 0  on 1 df,   p=1
## Wald test            = 0  on 1 df,   p=1
## Score (logrank) test = 0  on 1 df,   p=1
## 
## 
## 
## $DQ.Mismatch.No.
## $DQ.Mismatch.No.$Model
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##                     coef exp(coef) se(coef)      z     p
## DQ.Mismatch.No. -0.02535   0.97497  0.05212 -0.486 0.627
## 
## Likelihood ratio test=0.24  on 1 df, p=0.6235
## n= 98, number of events= 9 
## 
## $DQ.Mismatch.No.$Summary
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##                     coef exp(coef) se(coef)      z Pr(>|z|)
## DQ.Mismatch.No. -0.02535   0.97497  0.05212 -0.486    0.627
## 
##                 exp(coef) exp(-coef) lower .95 upper .95
## DQ.Mismatch.No.     0.975      1.026    0.8803      1.08
## 
## Concordance= 0.53  (se = 0.063 )
## Likelihood ratio test= 0.24  on 1 df,   p=0.6
## Wald test            = 0.24  on 1 df,   p=0.6
## Score (logrank) test = 0.24  on 1 df,   p=0.6
## 
## 
## 
## $age
## $age$Model
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##          coef exp(coef)  se(coef)      z     p
## age -0.002506  0.997497  0.029763 -0.084 0.933
## 
## Likelihood ratio test=0.01  on 1 df, p=0.9333
## n= 98, number of events= 9 
## 
## $age$Summary
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##          coef exp(coef)  se(coef)      z Pr(>|z|)
## age -0.002506  0.997497  0.029763 -0.084    0.933
## 
##     exp(coef) exp(-coef) lower .95 upper .95
## age    0.9975      1.003     0.941     1.057
## 
## Concordance= 0.477  (se = 0.111 )
## Likelihood ratio test= 0.01  on 1 df,   p=0.9
## Wald test            = 0.01  on 1 df,   p=0.9
## Score (logrank) test = 0.01  on 1 df,   p=0.9
cox_results2<-univariate_coxph(df, cumulative_cols, "status", "time2fail")

## 
## Cox proportional hazards model for ClassII....Mismatch..All :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##                              coef exp(coef) se(coef)     z Pr(>|z|)
## ClassII....Mismatch..All -0.01650   0.98364  0.01987 -0.83    0.406
## 
##                          exp(coef) exp(-coef) lower .95 upper .95
## ClassII....Mismatch..All    0.9836      1.017    0.9461     1.023
## 
## Concordance= 0.587  (se = 0.07 )
## Likelihood ratio test= 0.73  on 1 df,   p=0.4
## Wald test            = 0.69  on 1 df,   p=0.4
## Score (logrank) test = 0.68  on 1 df,   p=0.4

## 
## Cox proportional hazards model for ClassI....Mismatch..All :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##                            coef exp(coef) se(coef)     z Pr(>|z|)
## ClassI....Mismatch..All 0.03668   1.03736  0.02465 1.488    0.137
## 
##                         exp(coef) exp(-coef) lower .95 upper .95
## ClassI....Mismatch..All     1.037      0.964    0.9884     1.089
## 
## Concordance= 0.673  (se = 0.096 )
## Likelihood ratio test= 2.09  on 1 df,   p=0.1
## Wald test            = 2.21  on 1 df,   p=0.1
## Score (logrank) test = 2.2  on 1 df,   p=0.1

## 
## Cox proportional hazards model for tot_mismatch :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##                  coef exp(coef) se(coef)     z Pr(>|z|)
## tot_mismatch 0.002003  1.002005 0.014023 0.143    0.886
## 
##              exp(coef) exp(-coef) lower .95 upper .95
## tot_mismatch     1.002      0.998    0.9748      1.03
## 
## Concordance= 0.518  (se = 0.081 )
## Likelihood ratio test= 0.02  on 1 df,   p=0.9
## Wald test            = 0.02  on 1 df,   p=0.9
## Score (logrank) test = 0.02  on 1 df,   p=0.9
## 
##                                          Variable        HR  Lower_CI Upper_CI
## ClassII....Mismatch..All ClassII....Mismatch..All 0.9836383 0.9460761 1.022692
## ClassI....Mismatch..All   ClassI....Mismatch..All 1.0373620 0.9884388 1.088707
## tot_mismatch                         tot_mismatch 1.0020053 0.9748400 1.029927
##                            p_value   C_index
## ClassII....Mismatch..All 0.4062866 0.5873533
## ClassI....Mismatch..All  0.1367029 0.6734029
## tot_mismatch             0.8864073 0.5182529
cox_results2
## $ClassII....Mismatch..All
## $ClassII....Mismatch..All$Model
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##                              coef exp(coef) se(coef)     z     p
## ClassII....Mismatch..All -0.01650   0.98364  0.01987 -0.83 0.406
## 
## Likelihood ratio test=0.73  on 1 df, p=0.3939
## n= 98, number of events= 9 
## 
## $ClassII....Mismatch..All$Summary
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##                              coef exp(coef) se(coef)     z Pr(>|z|)
## ClassII....Mismatch..All -0.01650   0.98364  0.01987 -0.83    0.406
## 
##                          exp(coef) exp(-coef) lower .95 upper .95
## ClassII....Mismatch..All    0.9836      1.017    0.9461     1.023
## 
## Concordance= 0.587  (se = 0.07 )
## Likelihood ratio test= 0.73  on 1 df,   p=0.4
## Wald test            = 0.69  on 1 df,   p=0.4
## Score (logrank) test = 0.68  on 1 df,   p=0.4
## 
## 
## 
## $ClassI....Mismatch..All
## $ClassI....Mismatch..All$Model
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##                            coef exp(coef) se(coef)     z     p
## ClassI....Mismatch..All 0.03668   1.03736  0.02465 1.488 0.137
## 
## Likelihood ratio test=2.09  on 1 df, p=0.148
## n= 98, number of events= 9 
## 
## $ClassI....Mismatch..All$Summary
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##                            coef exp(coef) se(coef)     z Pr(>|z|)
## ClassI....Mismatch..All 0.03668   1.03736  0.02465 1.488    0.137
## 
##                         exp(coef) exp(-coef) lower .95 upper .95
## ClassI....Mismatch..All     1.037      0.964    0.9884     1.089
## 
## Concordance= 0.673  (se = 0.096 )
## Likelihood ratio test= 2.09  on 1 df,   p=0.1
## Wald test            = 2.21  on 1 df,   p=0.1
## Score (logrank) test = 2.2  on 1 df,   p=0.1
## 
## 
## 
## $tot_mismatch
## $tot_mismatch$Model
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##                  coef exp(coef) se(coef)     z     p
## tot_mismatch 0.002003  1.002005 0.014023 0.143 0.886
## 
## Likelihood ratio test=0.02  on 1 df, p=0.8866
## n= 98, number of events= 9 
## 
## $tot_mismatch$Summary
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##                  coef exp(coef) se(coef)     z Pr(>|z|)
## tot_mismatch 0.002003  1.002005 0.014023 0.143    0.886
## 
##              exp(coef) exp(-coef) lower .95 upper .95
## tot_mismatch     1.002      0.998    0.9748      1.03
## 
## Concordance= 0.518  (se = 0.081 )
## Likelihood ratio test= 0.02  on 1 df,   p=0.9
## Wald test            = 0.02  on 1 df,   p=0.9
## Score (logrank) test = 0.02  on 1 df,   p=0.9

Time to failure for antigen MM

cox_results<-univariate_coxph(df, independent_al_vars, "status", "time2fail")
## 
## Cox proportional hazards model for A :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##      coef exp(coef) se(coef)      z Pr(>|z|)
## A -0.3304    0.7186   0.6660 -0.496     0.62
## 
##   exp(coef) exp(-coef) lower .95 upper .95
## A    0.7186      1.392    0.1948     2.651
## 
## Concordance= 0.493  (se = 0.092 )
## Likelihood ratio test= 0.24  on 1 df,   p=0.6
## Wald test            = 0.25  on 1 df,   p=0.6
## Score (logrank) test = 0.25  on 1 df,   p=0.6
## Warning in coxph.fit(X, Y, istrat, offset, init, control, weights = weights, :
## Loglik converged before variable 1 ; coefficient may be infinite.

## 
## Cox proportional hazards model for B :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##        coef exp(coef)  se(coef)     z Pr(>|z|)
## B 1.936e+01 2.565e+08 1.064e+04 0.002    0.999
## 
##   exp(coef) exp(-coef) lower .95 upper .95
## B 256543565  3.898e-09         0       Inf
## 
## Concordance= 0.602  (se = 0.021 )
## Likelihood ratio test= 4.04  on 1 df,   p=0.04
## Wald test            = 0  on 1 df,   p=1
## Score (logrank) test = 2.27  on 1 df,   p=0.1

## 
## Cox proportional hazards model for C :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##     coef exp(coef) se(coef)     z Pr(>|z|)
## C 1.0649    2.9005   0.7584 1.404     0.16
## 
##   exp(coef) exp(-coef) lower .95 upper .95
## C       2.9     0.3448     0.656     12.82
## 
## Concordance= 0.646  (se = 0.062 )
## Likelihood ratio test= 2.54  on 1 df,   p=0.1
## Wald test            = 1.97  on 1 df,   p=0.2
## Score (logrank) test = 2.15  on 1 df,   p=0.1

## 
## Cox proportional hazards model for DRB1 :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##        coef exp(coef) se(coef)     z Pr(>|z|)
## DRB1 0.9040    2.4694   0.7626 1.185    0.236
## 
##      exp(coef) exp(-coef) lower .95 upper .95
## DRB1     2.469      0.405    0.5539     11.01
## 
## Concordance= 0.619  (se = 0.063 )
## Likelihood ratio test= 1.74  on 1 df,   p=0.2
## Wald test            = 1.41  on 1 df,   p=0.2
## Score (logrank) test = 1.5  on 1 df,   p=0.2

## 
## Cox proportional hazards model for DRB345 :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##            coef exp(coef) se(coef)      z Pr(>|z|)
## DRB345 -0.02847   0.97193  0.18092 -0.157    0.875
## 
##        exp(coef) exp(-coef) lower .95 upper .95
## DRB345    0.9719      1.029    0.6818     1.386
## 
## Concordance= 0.481  (se = 0.08 )
## Likelihood ratio test= 0.02  on 1 df,   p=0.9
## Wald test            = 0.02  on 1 df,   p=0.9
## Score (logrank) test = 0.02  on 1 df,   p=0.9

## 
## Cox proportional hazards model for DP :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##      coef exp(coef) se(coef)     z Pr(>|z|)  
## DP 0.5567    1.7449   0.3350 1.662   0.0966 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##    exp(coef) exp(-coef) lower .95 upper .95
## DP     1.745     0.5731    0.9049     3.365
## 
## Concordance= 0.614  (se = 0.094 )
## Likelihood ratio test= 3.03  on 1 df,   p=0.08
## Wald test            = 2.76  on 1 df,   p=0.1
## Score (logrank) test = 2.84  on 1 df,   p=0.09

## 
## Cox proportional hazards model for DQ :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##      coef exp(coef) se(coef)    z Pr(>|z|)
## DQ 0.2210    1.2473   0.3252 0.68    0.497
## 
##    exp(coef) exp(-coef) lower .95 upper .95
## DQ     1.247     0.8017    0.6594     2.359
## 
## Concordance= 0.561  (se = 0.081 )
## Likelihood ratio test= 0.48  on 1 df,   p=0.5
## Wald test            = 0.46  on 1 df,   p=0.5
## Score (logrank) test = 0.47  on 1 df,   p=0.5

## 
## Cox proportional hazards model for age :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 9 
## 
##          coef exp(coef)  se(coef)      z Pr(>|z|)
## age -0.002506  0.997497  0.029763 -0.084    0.933
## 
##     exp(coef) exp(-coef) lower .95 upper .95
## age    0.9975      1.003     0.941     1.057
## 
## Concordance= 0.477  (se = 0.111 )
## Likelihood ratio test= 0.01  on 1 df,   p=0.9
## Wald test            = 0.01  on 1 df,   p=0.9
## Score (logrank) test = 0.01  on 1 df,   p=0.9
## 
##        Variable           HR  Lower_CI  Upper_CI    p_value   C_index
## A             A 7.186378e-01 0.1947981  2.651156 0.61984442 0.4928292
## B             B 2.565436e+08 0.0000000       Inf 0.99854768 0.6016949
## C             C 2.900495e+00 0.6559927 12.824645 0.16029853 0.6460235
## DRB1       DRB1 2.469357e+00 0.5539330 11.008053 0.23587424 0.6186441
## DRB345   DRB345 9.719275e-01 0.6817626  1.385590 0.87494277 0.4810952
## DP           DP 1.744910e+00 0.9049047  3.364676 0.09657353 0.6140808
## DQ           DQ 1.247328e+00 0.6594318  2.359344 0.49676253 0.5612777
## age         age 9.974973e-01 0.9409729  1.057417 0.93290381 0.4771838

Grade 1 coxph eplet

cox_results<-univariate_coxph(df, independent_vars, "Status.G1", "Med.Time.G1")
## Warning in splines::ns(temp, df = df, intercept = TRUE): shoving 'interior'
## knots matching boundary knots to inside
## 
## Cox proportional hazards model for A.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 96 
## 
##                   coef exp(coef) se(coef)     z Pr(>|z|)
## A.Mismatch.No. 0.00244   1.00244  0.01379 0.177     0.86
## 
##                exp(coef) exp(-coef) lower .95 upper .95
## A.Mismatch.No.     1.002     0.9976    0.9757      1.03
## 
## Concordance= 0.527  (se = 0.047 )
## Likelihood ratio test= 0.03  on 1 df,   p=0.9
## Wald test            = 0.03  on 1 df,   p=0.9
## Score (logrank) test = 0.03  on 1 df,   p=0.9
## Warning in splines::ns(temp, df = df, intercept = TRUE): shoving 'interior'
## knots matching boundary knots to inside

## 
## Cox proportional hazards model for B.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 96 
## 
##                   coef exp(coef) se(coef)     z Pr(>|z|)
## B.Mismatch.No. 0.02007   1.02028  0.01742 1.152    0.249
## 
##                exp(coef) exp(-coef) lower .95 upper .95
## B.Mismatch.No.      1.02     0.9801     0.986     1.056
## 
## Concordance= 0.547  (se = 0.044 )
## Likelihood ratio test= 1.32  on 1 df,   p=0.3
## Wald test            = 1.33  on 1 df,   p=0.2
## Score (logrank) test = 1.33  on 1 df,   p=0.2
## Warning in splines::ns(temp, df = df, intercept = TRUE): shoving 'interior'
## knots matching boundary knots to inside

## 
## Cox proportional hazards model for C.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 96 
## 
##                    coef exp(coef) se(coef)      z Pr(>|z|)
## C.Mismatch.No. -0.01988   0.98032  0.02019 -0.985    0.325
## 
##                exp(coef) exp(-coef) lower .95 upper .95
## C.Mismatch.No.    0.9803       1.02    0.9423      1.02
## 
## Concordance= 0.56  (se = 0.047 )
## Likelihood ratio test= 0.98  on 1 df,   p=0.3
## Wald test            = 0.97  on 1 df,   p=0.3
## Score (logrank) test = 0.97  on 1 df,   p=0.3
## Warning in splines::ns(temp, df = df, intercept = TRUE): shoving 'interior'
## knots matching boundary knots to inside

## 
## Cox proportional hazards model for DRB1.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 96 
## 
##                       coef exp(coef) se(coef)     z Pr(>|z|)
## DRB1.Mismatch.No. -0.01522   0.98490  0.01537 -0.99    0.322
## 
##                   exp(coef) exp(-coef) lower .95 upper .95
## DRB1.Mismatch.No.    0.9849      1.015    0.9557     1.015
## 
## Concordance= 0.506  (se = 0.047 )
## Likelihood ratio test= 0.98  on 1 df,   p=0.3
## Wald test            = 0.98  on 1 df,   p=0.3
## Score (logrank) test = 0.98  on 1 df,   p=0.3
## Warning in splines::ns(temp, df = df, intercept = TRUE): shoving 'interior'
## knots matching boundary knots to inside

## 
## Cox proportional hazards model for DRB345.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 96 
## 
##                         coef exp(coef) se(coef)      z Pr(>|z|)  
## DRB345.Mismatch.No. -0.03434   0.96625  0.01359 -2.527   0.0115 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##                     exp(coef) exp(-coef) lower .95 upper .95
## DRB345.Mismatch.No.    0.9662      1.035    0.9409    0.9923
## 
## Concordance= 0.602  (se = 0.045 )
## Likelihood ratio test= 6.74  on 1 df,   p=0.009
## Wald test            = 6.39  on 1 df,   p=0.01
## Score (logrank) test = 6.36  on 1 df,   p=0.01
## Warning in splines::ns(temp, df = df, intercept = TRUE): shoving 'interior'
## knots matching boundary knots to inside

## 
## Cox proportional hazards model for DP.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 96 
## 
##                      coef exp(coef)  se(coef)     z Pr(>|z|)
## DP.Mismatch.No. -0.006597  0.993425  0.014647 -0.45    0.652
## 
##                 exp(coef) exp(-coef) lower .95 upper .95
## DP.Mismatch.No.    0.9934      1.007    0.9653     1.022
## 
## Concordance= 0.545  (se = 0.049 )
## Likelihood ratio test= 0.2  on 1 df,   p=0.7
## Wald test            = 0.2  on 1 df,   p=0.7
## Score (logrank) test = 0.2  on 1 df,   p=0.7
## Warning in splines::ns(temp, df = df, intercept = TRUE): shoving 'interior'
## knots matching boundary knots to inside

## 
## Cox proportional hazards model for DQ.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 96 
## 
##                    coef exp(coef) se(coef)     z Pr(>|z|)
## DQ.Mismatch.No. 0.01773   1.01789  0.01536 1.155    0.248
## 
##                 exp(coef) exp(-coef) lower .95 upper .95
## DQ.Mismatch.No.     1.018     0.9824    0.9877     1.049
## 
## Concordance= 0.585  (se = 0.041 )
## Likelihood ratio test= 1.31  on 1 df,   p=0.3
## Wald test            = 1.33  on 1 df,   p=0.2
## Score (logrank) test = 1.33  on 1 df,   p=0.2
## Warning in splines::ns(temp, df = df, intercept = TRUE): shoving 'interior'
## knots matching boundary knots to inside

## 
## Cox proportional hazards model for age :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 96 
## 
##         coef exp(coef) se(coef)      z Pr(>|z|)    
## age -0.03635   0.96430  0.00961 -3.783 0.000155 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##     exp(coef) exp(-coef) lower .95 upper .95
## age    0.9643      1.037    0.9463    0.9826
## 
## Concordance= 0.638  (se = 0.038 )
## Likelihood ratio test= 12.55  on 1 df,   p=4e-04
## Wald test            = 14.31  on 1 df,   p=2e-04
## Score (logrank) test = 14.54  on 1 df,   p=1e-04
## 
##                                Variable        HR  Lower_CI Upper_CI
## A.Mismatch.No.           A.Mismatch.No. 1.0024427 0.9757186 1.029899
## B.Mismatch.No.           B.Mismatch.No. 1.0202758 0.9860249 1.055717
## C.Mismatch.No.           C.Mismatch.No. 0.9803180 0.9422908 1.019880
## DRB1.Mismatch.No.     DRB1.Mismatch.No. 0.9848956 0.9556609 1.015025
## DRB345.Mismatch.No. DRB345.Mismatch.No. 0.9662456 0.9408553 0.992321
## DP.Mismatch.No.         DP.Mismatch.No. 0.9934252 0.9653119 1.022357
## DQ.Mismatch.No.         DQ.Mismatch.No. 1.0178887 0.9877095 1.048990
## age                                 age 0.9643023 0.9463096 0.982637
##                          p_value   C_index
## A.Mismatch.No.      0.8595337966 0.5266052
## B.Mismatch.No.      0.2492568889 0.5466706
## C.Mismatch.No.      0.3247370730 0.5600476
## DRB1.Mismatch.No.   0.3221961439 0.5059453
## DRB345.Mismatch.No. 0.0114928267 0.6015161
## DP.Mismatch.No.     0.6524459200 0.5447384
## DQ.Mismatch.No.     0.2482394454 0.5847206
## age                 0.0001551797 0.6382283
cox_multi<-coxph(Surv(Med.Time.G1, Status.G1) ~ A.Mismatch.No. + B.Mismatch.No. + C.Mismatch.No. + DRB1.Mismatch.No. + DRB345.Mismatch.No. + DP.Mismatch.No. + DQ.Mismatch.No. + age, data=df)
summary(cox_multi)
## Call:
## coxph(formula = Surv(Med.Time.G1, Status.G1) ~ A.Mismatch.No. + 
##     B.Mismatch.No. + C.Mismatch.No. + DRB1.Mismatch.No. + DRB345.Mismatch.No. + 
##     DP.Mismatch.No. + DQ.Mismatch.No. + age, data = df)
## 
##   n= 98, number of events= 96 
## 
##                          coef exp(coef)  se(coef)      z Pr(>|z|)   
## A.Mismatch.No.       0.005557  1.005573  0.016575  0.335  0.73742   
## B.Mismatch.No.       0.015307  1.015425  0.019590  0.781  0.43460   
## C.Mismatch.No.      -0.034148  0.966428  0.024514 -1.393  0.16362   
## DRB1.Mismatch.No.    0.007341  1.007368  0.020268  0.362  0.71720   
## DRB345.Mismatch.No. -0.043407  0.957522  0.016613 -2.613  0.00898 **
## DP.Mismatch.No.     -0.018227  0.981938  0.015131 -1.205  0.22835   
## DQ.Mismatch.No.      0.031224  1.031716  0.018178  1.718  0.08586 . 
## age                 -0.031778  0.968722  0.010238 -3.104  0.00191 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##                     exp(coef) exp(-coef) lower .95 upper .95
## A.Mismatch.No.         1.0056     0.9945    0.9734    1.0388
## B.Mismatch.No.         1.0154     0.9848    0.9772    1.0552
## C.Mismatch.No.         0.9664     1.0347    0.9211    1.0140
## DRB1.Mismatch.No.      1.0074     0.9927    0.9681    1.0482
## DRB345.Mismatch.No.    0.9575     1.0444    0.9268    0.9892
## DP.Mismatch.No.        0.9819     1.0184    0.9532    1.0115
## DQ.Mismatch.No.        1.0317     0.9693    0.9956    1.0691
## age                    0.9687     1.0323    0.9495    0.9884
## 
## Concordance= 0.727  (se = 0.035 )
## Likelihood ratio test= 24.79  on 8 df,   p=0.002
## Wald test            = 25.71  on 8 df,   p=0.001
## Score (logrank) test = 26.9  on 8 df,   p=7e-04

Grade 1 coxph antigen

cox_results<-univariate_coxph(df, independent_al_vars, "Status.G1", "Med.Time.G1")
## Warning in splines::ns(temp, df = df, intercept = TRUE): shoving 'interior'
## knots matching boundary knots to inside
## 
## Cox proportional hazards model for A :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 96 
## 
##     coef exp(coef) se(coef)     z Pr(>|z|)
## A 0.1665    1.1811   0.1899 0.876    0.381
## 
##   exp(coef) exp(-coef) lower .95 upper .95
## A     1.181     0.8467     0.814     1.714
## 
## Concordance= 0.516  (se = 0.039 )
## Likelihood ratio test= 0.79  on 1 df,   p=0.4
## Wald test            = 0.77  on 1 df,   p=0.4
## Score (logrank) test = 0.77  on 1 df,   p=0.4
## Warning in splines::ns(temp, df = df, intercept = TRUE): shoving 'interior'
## knots matching boundary knots to inside

## 
## Cox proportional hazards model for B :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 96 
## 
##      coef exp(coef) se(coef)     z Pr(>|z|)
## B 0.06843   1.07083  0.25734 0.266     0.79
## 
##   exp(coef) exp(-coef) lower .95 upper .95
## B     1.071     0.9339    0.6466     1.773
## 
## Concordance= 0.521  (se = 0.031 )
## Likelihood ratio test= 0.07  on 1 df,   p=0.8
## Wald test            = 0.07  on 1 df,   p=0.8
## Score (logrank) test = 0.07  on 1 df,   p=0.8
## Warning in splines::ns(temp, df = df, intercept = TRUE): shoving 'interior'
## knots matching boundary knots to inside

## 
## Cox proportional hazards model for C :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 96 
## 
##      coef exp(coef) se(coef)      z Pr(>|z|)
## C -0.1132    0.8930   0.1589 -0.713    0.476
## 
##   exp(coef) exp(-coef) lower .95 upper .95
## C     0.893       1.12    0.6541     1.219
## 
## Concordance= 0.528  (se = 0.041 )
## Likelihood ratio test= 0.5  on 1 df,   p=0.5
## Wald test            = 0.51  on 1 df,   p=0.5
## Score (logrank) test = 0.51  on 1 df,   p=0.5
## Warning in splines::ns(temp, df = df, intercept = TRUE): shoving 'interior'
## knots matching boundary knots to inside

## 
## Cox proportional hazards model for DRB1 :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 96 
## 
##          coef exp(coef) se(coef)      z Pr(>|z|)
## DRB1 -0.00464   0.99537  0.17558 -0.026    0.979
## 
##      exp(coef) exp(-coef) lower .95 upper .95
## DRB1    0.9954      1.005    0.7056     1.404
## 
## Concordance= 0.483  (se = 0.038 )
## Likelihood ratio test= 0  on 1 df,   p=1
## Wald test            = 0  on 1 df,   p=1
## Score (logrank) test = 0  on 1 df,   p=1
## Warning in splines::ns(temp, df = df, intercept = TRUE): shoving 'interior'
## knots matching boundary knots to inside

## 
## Cox proportional hazards model for DRB345 :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 96 
## 
##            coef exp(coef) se(coef)      z Pr(>|z|)
## DRB345 -0.05932   0.94241  0.05821 -1.019    0.308
## 
##        exp(coef) exp(-coef) lower .95 upper .95
## DRB345    0.9424      1.061    0.8408     1.056
## 
## Concordance= 0.518  (se = 0.042 )
## Likelihood ratio test= 1.05  on 1 df,   p=0.3
## Wald test            = 1.04  on 1 df,   p=0.3
## Score (logrank) test = 1.04  on 1 df,   p=0.3
## Warning in splines::ns(temp, df = df, intercept = TRUE): shoving 'interior'
## knots matching boundary knots to inside

## 
## Cox proportional hazards model for DP :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 96 
## 
##        coef exp(coef) se(coef)      z Pr(>|z|)  
## DP -0.15745   0.85432  0.09407 -1.674   0.0942 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##    exp(coef) exp(-coef) lower .95 upper .95
## DP    0.8543      1.171    0.7105     1.027
## 
## Concordance= 0.582  (se = 0.045 )
## Likelihood ratio test= 2.78  on 1 df,   p=0.1
## Wald test            = 2.8  on 1 df,   p=0.09
## Score (logrank) test = 2.81  on 1 df,   p=0.09
## Warning in splines::ns(temp, df = df, intercept = TRUE): shoving 'interior'
## knots matching boundary knots to inside

## 
## Cox proportional hazards model for DQ :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 96 
## 
##       coef exp(coef) se(coef)   z Pr(>|z|)
## DQ 0.04681   1.04793  0.09354 0.5    0.617
## 
##    exp(coef) exp(-coef) lower .95 upper .95
## DQ     1.048     0.9543    0.8724     1.259
## 
## Concordance= 0.547  (se = 0.045 )
## Likelihood ratio test= 0.25  on 1 df,   p=0.6
## Wald test            = 0.25  on 1 df,   p=0.6
## Score (logrank) test = 0.25  on 1 df,   p=0.6
## Warning in splines::ns(temp, df = df, intercept = TRUE): shoving 'interior'
## knots matching boundary knots to inside

## 
## Cox proportional hazards model for age :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 96 
## 
##         coef exp(coef) se(coef)      z Pr(>|z|)    
## age -0.03635   0.96430  0.00961 -3.783 0.000155 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##     exp(coef) exp(-coef) lower .95 upper .95
## age    0.9643      1.037    0.9463    0.9826
## 
## Concordance= 0.638  (se = 0.038 )
## Likelihood ratio test= 12.55  on 1 df,   p=4e-04
## Wald test            = 14.31  on 1 df,   p=2e-04
## Score (logrank) test = 14.54  on 1 df,   p=1e-04
## 
##        Variable        HR  Lower_CI Upper_CI      p_value   C_index
## A             A 1.1811174 0.8140199 1.713764 0.3807632086 0.5157551
## B             B 1.0708253 0.6466464 1.773252 0.7903107441 0.5212545
## C             C 0.8929849 0.6540734 1.219163 0.4761499600 0.5280916
## DRB1       DRB1 0.9953712 0.7055544 1.404234 0.9789192987 0.4829073
## DRB345   DRB345 0.9424080 0.8407948 1.056302 0.3081997742 0.5184304
## DP           DP 0.8543197 0.7104707 1.027294 0.0941863894 0.5820452
## DQ           DQ 1.0479271 0.8723862 1.258790 0.6167478034 0.5466706
## age         age 0.9643023 0.9463096 0.982637 0.0001551797 0.6382283
cox_multi<-coxph(Surv(Med.Time.G1, Status.G1) ~ A + B + C + DRB1 + DRB345 + DP + DQ +age, data=df)
summary(cox_multi)
## Call:
## coxph(formula = Surv(Med.Time.G1, Status.G1) ~ A + B + C + DRB1 + 
##     DRB345 + DP + DQ + age, data = df)
## 
##   n= 98, number of events= 96 
## 
##            coef exp(coef) se(coef)      z Pr(>|z|)    
## A       0.15483   1.16746  0.20462  0.757  0.44926    
## B       0.15349   1.16590  0.28542  0.538  0.59073    
## C      -0.12841   0.87949  0.17697 -0.726  0.46808    
## DRB1    0.01155   1.01162  0.23437  0.049  0.96068    
## DRB345 -0.06878   0.93353  0.06894 -0.998  0.31846    
## DP     -0.19582   0.82216  0.09935 -1.971  0.04872 *  
## DQ      0.08192   1.08536  0.11669  0.702  0.48268    
## age    -0.03373   0.96683  0.01001 -3.371  0.00075 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##        exp(coef) exp(-coef) lower .95 upper .95
## A         1.1675     0.8566    0.7817    1.7435
## B         1.1659     0.8577    0.6664    2.0399
## C         0.8795     1.1370    0.6217    1.2441
## DRB1      1.0116     0.9885    0.6390    1.6014
## DRB345    0.9335     1.0712    0.8155    1.0686
## DP        0.8222     1.2163    0.6767    0.9989
## DQ        1.0854     0.9213    0.8635    1.3643
## age       0.9668     1.0343    0.9481    0.9860
## 
## Concordance= 0.696  (se = 0.037 )
## Likelihood ratio test= 18.34  on 8 df,   p=0.02
## Wald test            = 20.05  on 8 df,   p=0.01
## Score (logrank) test = 20.55  on 8 df,   p=0.008
ggcoxzph(cox.zph(cox_multi))
## Warning in splines::ns(temp, df = df, intercept = TRUE): shoving 'interior'
## knots matching boundary knots to inside

Grade 2 coxph eplet

cox_results<-univariate_coxph(df, independent_vars, "Status.G2", "Med.Time.G2")

## 
## Cox proportional hazards model for A.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 34 
## 
##                   coef exp(coef) se(coef)     z Pr(>|z|)
## A.Mismatch.No. 0.03169   1.03220  0.02436 1.301    0.193
## 
##                exp(coef) exp(-coef) lower .95 upper .95
## A.Mismatch.No.     1.032     0.9688    0.9841     1.083
## 
## Concordance= 0.581  (se = 0.051 )
## Likelihood ratio test= 1.69  on 1 df,   p=0.2
## Wald test            = 1.69  on 1 df,   p=0.2
## Score (logrank) test = 1.69  on 1 df,   p=0.2

## 
## Cox proportional hazards model for B.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 34 
## 
##                   coef exp(coef) se(coef)     z Pr(>|z|)  
## B.Mismatch.No. 0.04873   1.04994  0.02911 1.674   0.0941 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##                exp(coef) exp(-coef) lower .95 upper .95
## B.Mismatch.No.      1.05     0.9524    0.9917     1.112
## 
## Concordance= 0.581  (se = 0.053 )
## Likelihood ratio test= 2.8  on 1 df,   p=0.09
## Wald test            = 2.8  on 1 df,   p=0.09
## Score (logrank) test = 2.82  on 1 df,   p=0.09

## 
## Cox proportional hazards model for C.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 34 
## 
##                    coef exp(coef) se(coef)     z Pr(>|z|)
## C.Mismatch.No. 0.006271  1.006290 0.033584 0.187    0.852
## 
##                exp(coef) exp(-coef) lower .95 upper .95
## C.Mismatch.No.     1.006     0.9937    0.9422     1.075
## 
## Concordance= 0.508  (se = 0.051 )
## Likelihood ratio test= 0.03  on 1 df,   p=0.9
## Wald test            = 0.03  on 1 df,   p=0.9
## Score (logrank) test = 0.03  on 1 df,   p=0.9

## 
## Cox proportional hazards model for DRB1.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 34 
## 
##                      coef exp(coef) se(coef)     z Pr(>|z|)
## DRB1.Mismatch.No. 0.04261   1.04353  0.02795 1.525    0.127
## 
##                   exp(coef) exp(-coef) lower .95 upper .95
## DRB1.Mismatch.No.     1.044     0.9583    0.9879     1.102
## 
## Concordance= 0.577  (se = 0.054 )
## Likelihood ratio test= 2.32  on 1 df,   p=0.1
## Wald test            = 2.33  on 1 df,   p=0.1
## Score (logrank) test = 2.34  on 1 df,   p=0.1

## 
## Cox proportional hazards model for DRB345.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 34 
## 
##                         coef exp(coef) se(coef)     z Pr(>|z|)
## DRB345.Mismatch.No. 0.006783  1.006806 0.020400 0.333    0.739
## 
##                     exp(coef) exp(-coef) lower .95 upper .95
## DRB345.Mismatch.No.     1.007     0.9932    0.9673     1.048
## 
## Concordance= 0.537  (se = 0.056 )
## Likelihood ratio test= 0.11  on 1 df,   p=0.7
## Wald test            = 0.11  on 1 df,   p=0.7
## Score (logrank) test = 0.11  on 1 df,   p=0.7

## 
## Cox proportional hazards model for DP.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 34 
## 
##                     coef exp(coef) se(coef)      z Pr(>|z|)
## DP.Mismatch.No. -0.03786   0.96285  0.02566 -1.475     0.14
## 
##                 exp(coef) exp(-coef) lower .95 upper .95
## DP.Mismatch.No.    0.9628      1.039    0.9156     1.013
## 
## Concordance= 0.584  (se = 0.054 )
## Likelihood ratio test= 2.26  on 1 df,   p=0.1
## Wald test            = 2.18  on 1 df,   p=0.1
## Score (logrank) test = 2.2  on 1 df,   p=0.1

## 
## Cox proportional hazards model for DQ.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 34 
## 
##                    coef exp(coef) se(coef)     z Pr(>|z|)   
## DQ.Mismatch.No. 0.06890   1.07133  0.02399 2.872  0.00408 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##                 exp(coef) exp(-coef) lower .95 upper .95
## DQ.Mismatch.No.     1.071     0.9334     1.022     1.123
## 
## Concordance= 0.661  (se = 0.046 )
## Likelihood ratio test= 7.87  on 1 df,   p=0.005
## Wald test            = 8.25  on 1 df,   p=0.004
## Score (logrank) test = 8.36  on 1 df,   p=0.004

## 
## Cox proportional hazards model for age :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 34 
## 
##         coef exp(coef) se(coef)      z Pr(>|z|)   
## age -0.04138   0.95946  0.01317 -3.141  0.00168 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##     exp(coef) exp(-coef) lower .95 upper .95
## age    0.9595      1.042     0.935    0.9846
## 
## Concordance= 0.631  (se = 0.054 )
## Likelihood ratio test= 8.29  on 1 df,   p=0.004
## Wald test            = 9.87  on 1 df,   p=0.002
## Score (logrank) test = 10.14  on 1 df,   p=0.001
## 
##                                Variable        HR  Lower_CI  Upper_CI
## A.Mismatch.No.           A.Mismatch.No. 1.0322016 0.9840738 1.0826832
## B.Mismatch.No.           B.Mismatch.No. 1.0499368 0.9917195 1.1115716
## C.Mismatch.No.           C.Mismatch.No. 1.0062902 0.9421854 1.0747567
## DRB1.Mismatch.No.     DRB1.Mismatch.No. 1.0435337 0.9879151 1.1022835
## DRB345.Mismatch.No. DRB345.Mismatch.No. 1.0068064 0.9673459 1.0478767
## DP.Mismatch.No.         DP.Mismatch.No. 0.9628452 0.9156111 1.0125160
## DQ.Mismatch.No.         DQ.Mismatch.No. 1.0713328 1.0221187 1.1229166
## age                                 age 0.9594634 0.9350077 0.9845588
##                         p_value   C_index
## A.Mismatch.No.      0.193270122 0.5811124
## B.Mismatch.No.      0.094076202 0.5813055
## C.Mismatch.No.      0.851887464 0.5075319
## DRB1.Mismatch.No.   0.127289682 0.5770568
## DRB345.Mismatch.No. 0.739492864 0.5365006
## DP.Mismatch.No.     0.140129484 0.5836230
## DQ.Mismatch.No.     0.004081578 0.6612592
## age                 0.001682256 0.6307455
cox_multi<-coxph(Surv(Med.Time.G2, Status.G2) ~ A.Mismatch.No. + B.Mismatch.No. + C.Mismatch.No. + DRB1.Mismatch.No. + DRB345.Mismatch.No. + DP.Mismatch.No. + DQ.Mismatch.No. +age, data=df)
ggcoxzph(cox.zph(cox_multi))

summary(cox_multi)
## Call:
## coxph(formula = Surv(Med.Time.G2, Status.G2) ~ A.Mismatch.No. + 
##     B.Mismatch.No. + C.Mismatch.No. + DRB1.Mismatch.No. + DRB345.Mismatch.No. + 
##     DP.Mismatch.No. + DQ.Mismatch.No. + age, data = df)
## 
##   n= 98, number of events= 34 
## 
##                          coef exp(coef)  se(coef)      z Pr(>|z|)    
## A.Mismatch.No.       0.008364  1.008400  0.027652  0.302 0.762276    
## B.Mismatch.No.       0.024613  1.024918  0.032166  0.765 0.444161    
## C.Mismatch.No.       0.003002  1.003006  0.039270  0.076 0.939075    
## DRB1.Mismatch.No.    0.047584  1.048734  0.036248  1.313 0.189275    
## DRB345.Mismatch.No. -0.012720  0.987360  0.023776 -0.535 0.592641    
## DP.Mismatch.No.     -0.062483  0.939429  0.025554 -2.445 0.014478 *  
## DQ.Mismatch.No.      0.058841  1.060606  0.026983  2.181 0.029212 *  
## age                 -0.054905  0.946575  0.015141 -3.626 0.000287 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##                     exp(coef) exp(-coef) lower .95 upper .95
## A.Mismatch.No.         1.0084     0.9917    0.9552    1.0646
## B.Mismatch.No.         1.0249     0.9757    0.9623    1.0916
## C.Mismatch.No.         1.0030     0.9970    0.9287    1.0833
## DRB1.Mismatch.No.      1.0487     0.9535    0.9768    1.1260
## DRB345.Mismatch.No.    0.9874     1.0128    0.9424    1.0345
## DP.Mismatch.No.        0.9394     1.0645    0.8935    0.9877
## DQ.Mismatch.No.        1.0606     0.9429    1.0060    1.1182
## age                    0.9466     1.0564    0.9189    0.9751
## 
## Concordance= 0.744  (se = 0.046 )
## Likelihood ratio test= 24.79  on 8 df,   p=0.002
## Wald test            = 24.82  on 8 df,   p=0.002
## Score (logrank) test = 25.7  on 8 df,   p=0.001

Grade 2 coxph antigen

cox_results<-univariate_coxph(df, independent_al_vars, "Status.G2", "Med.Time.G2")

## 
## Cox proportional hazards model for A :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 34 
## 
##      coef exp(coef) se(coef)     z Pr(>|z|)
## A 0.09735   1.10225  0.34331 0.284    0.777
## 
##   exp(coef) exp(-coef) lower .95 upper .95
## A     1.102     0.9072    0.5624      2.16
## 
## Concordance= 0.506  (se = 0.044 )
## Likelihood ratio test= 0.08  on 1 df,   p=0.8
## Wald test            = 0.08  on 1 df,   p=0.8
## Score (logrank) test = 0.08  on 1 df,   p=0.8

## 
## Cox proportional hazards model for B :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 34 
## 
##    coef exp(coef) se(coef)     z Pr(>|z|)  
## B 1.096     2.993    0.605 1.812     0.07 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##   exp(coef) exp(-coef) lower .95 upper .95
## B     2.992     0.3342    0.9143     9.795
## 
## Concordance= 0.575  (se = 0.029 )
## Likelihood ratio test= 4.43  on 1 df,   p=0.04
## Wald test            = 3.28  on 1 df,   p=0.07
## Score (logrank) test = 3.62  on 1 df,   p=0.06

## 
## Cox proportional hazards model for C :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 34 
## 
##      coef exp(coef) se(coef)     z Pr(>|z|)
## C 0.08003   1.08332  0.28292 0.283    0.777
## 
##   exp(coef) exp(-coef) lower .95 upper .95
## C     1.083     0.9231    0.6222     1.886
## 
## Concordance= 0.514  (se = 0.045 )
## Likelihood ratio test= 0.08  on 1 df,   p=0.8
## Wald test            = 0.08  on 1 df,   p=0.8
## Score (logrank) test = 0.08  on 1 df,   p=0.8

## 
## Cox proportional hazards model for DRB1 :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 34 
## 
##        coef exp(coef) se(coef)     z Pr(>|z|)
## DRB1 0.3188    1.3754   0.3196 0.997    0.319
## 
##      exp(coef) exp(-coef) lower .95 upper .95
## DRB1     1.375      0.727    0.7351     2.573
## 
## Concordance= 0.546  (se = 0.043 )
## Likelihood ratio test= 1.06  on 1 df,   p=0.3
## Wald test            = 0.99  on 1 df,   p=0.3
## Score (logrank) test = 1  on 1 df,   p=0.3

## 
## Cox proportional hazards model for DRB345 :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 34 
## 
##             coef exp(coef)  se(coef)      z Pr(>|z|)
## DRB345 -0.005122  0.994891  0.093669 -0.055    0.956
## 
##        exp(coef) exp(-coef) lower .95 upper .95
## DRB345    0.9949      1.005     0.828     1.195
## 
## Concordance= 0.491  (se = 0.05 )
## Likelihood ratio test= 0  on 1 df,   p=1
## Wald test            = 0  on 1 df,   p=1
## Score (logrank) test = 0  on 1 df,   p=1

## 
## Cox proportional hazards model for DP :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 34 
## 
##       coef exp(coef) se(coef)      z Pr(>|z|)  
## DP -0.4018    0.6691   0.1568 -2.562   0.0104 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##    exp(coef) exp(-coef) lower .95 upper .95
## DP    0.6691      1.494    0.4921    0.9099
## 
## Concordance= 0.63  (se = 0.053 )
## Likelihood ratio test= 6.48  on 1 df,   p=0.01
## Wald test            = 6.57  on 1 df,   p=0.01
## Score (logrank) test = 6.72  on 1 df,   p=0.01

## 
## Cox proportional hazards model for DQ :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 34 
## 
##      coef exp(coef) se(coef)     z Pr(>|z|)  
## DQ 0.4520    1.5715   0.1809 2.499   0.0124 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##    exp(coef) exp(-coef) lower .95 upper .95
## DQ     1.571     0.6363     1.102      2.24
## 
## Concordance= 0.624  (se = 0.049 )
## Likelihood ratio test= 6.84  on 1 df,   p=0.009
## Wald test            = 6.25  on 1 df,   p=0.01
## Score (logrank) test = 6.41  on 1 df,   p=0.01

## 
## Cox proportional hazards model for age :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 98, number of events= 34 
## 
##         coef exp(coef) se(coef)      z Pr(>|z|)   
## age -0.04138   0.95946  0.01317 -3.141  0.00168 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##     exp(coef) exp(-coef) lower .95 upper .95
## age    0.9595      1.042     0.935    0.9846
## 
## Concordance= 0.631  (se = 0.054 )
## Likelihood ratio test= 8.29  on 1 df,   p=0.004
## Wald test            = 9.87  on 1 df,   p=0.002
## Score (logrank) test = 10.14  on 1 df,   p=0.001
## 
##        Variable        HR  Lower_CI  Upper_CI     p_value   C_index
## A             A 1.1022508 0.5624173 2.1602409 0.776732635 0.5057937
## B             B 2.9924740 0.9142523 9.7947807 0.070020100 0.5749324
## C             C 1.0833205 0.6222100 1.8861531 0.777271200 0.5142912
## DRB1       DRB1 1.3754244 0.7351120 2.5734750 0.318650296 0.5459637
## DRB345   DRB345 0.9948913 0.8280269 1.1953822 0.956393252 0.4911163
## DP           DP 0.6691286 0.4920872 0.9098652 0.010395637 0.6295867
## DQ           DQ 1.5714794 1.1024534 2.2400471 0.012445126 0.6235998
## age         age 0.9594634 0.9350077 0.9845588 0.001682256 0.6307455

data pivoting longer

time_cols<-c("X1", "X2", "X3", "X4", "X6", "X8", "X10", "X12", "X26", "X52")
df_long<-df %>% pivot_longer(cols = time_cols, names_to = "Time_biopsy", values_to = "Grade")
## Warning: Using an external vector in selections was deprecated in tidyselect 1.1.0.
## ℹ Please use `all_of()` or `any_of()` instead.
##   # Was:
##   data %>% select(time_cols)
## 
##   # Now:
##   data %>% select(all_of(time_cols))
## 
## See <https://tidyselect.r-lib.org/reference/faq-external-vector.html>.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
df_long$Time_biopsy<-df_long$Time_biopsy %>% gsub("X", "", .) %>% as.numeric()
df_long$Grade[df_long$Grade==9]<-NA
df_long$Grade<-df_long$Grade+1
df_long$Grade<-df_long$Grade %>% as.factor()
df_long <- df_long[order(df_long$Time_biopsy),]
df_long<-df_long %>% filter(!is.na(Grade))

function to write repeated measures PO model

create_pomtim_formula <- function(response_var, independent_vars, time,df) {
  all_vars <- c(response_var, independent_vars) 
  missing_vars <- setdiff(all_vars, names(df))
  
  if(length(missing_vars) > 0) { stop(paste("The following variables are missing in the dataframe:",
                                             paste(missing_vars, collapse=", "), ".")) 
    }
  formula_str <- paste(response_var, "~", paste(independent_vars, collapse = " + "), "+", time)
  return(as.formula(formula_str))
}
create_pom_formula <- function(response_var, independent_vars, time,df) {
  all_vars <- c(response_var, independent_vars) 
  missing_vars <- setdiff(all_vars, names(df))
  
  if(length(missing_vars) > 0) { stop(paste("The following variables are missing in the dataframe:",
                                             paste(missing_vars, collapse=", "), ".")) 
    }
  formula_str <- paste(response_var, "~", paste(independent_vars, collapse = " + "))
  return(as.formula(formula_str))
}

Repeated measures PO model

independent_vars
## [1] "A.Mismatch.No."      "B.Mismatch.No."      "C.Mismatch.No."     
## [4] "DRB1.Mismatch.No."   "DRB345.Mismatch.No." "DP.Mismatch.No."    
## [7] "DQ.Mismatch.No."     "age"
pom_formula <- create_pomtim_formula("Grade", independent_vars, "Time_biopsy" ,df_long)
pom_formula
## Grade ~ A.Mismatch.No. + B.Mismatch.No. + C.Mismatch.No. + DRB1.Mismatch.No. + 
##     DRB345.Mismatch.No. + DP.Mismatch.No. + DQ.Mismatch.No. + 
##     age + Time_biopsy
## <environment: 0x5f649afbaa80>
model <- polr(pom_formula, data = df_long)
#model <- polr(Grade ~ Time_biopsy, data = df_long, Hess = TRUE)
coefficients <- summary(model)$coefficients
## 
## Re-fitting to get Hessian
p_value <- (1 - pnorm(abs(coefficients[ ,"t value"]), 0, 1))*2

coefficients <- summary(model)$coefficients
## 
## Re-fitting to get Hessian
coefficients <- cbind(coefficients, p_value)

# calculate odds ratios
coefficients <- cbind(coefficients, odds_ratio = exp(coefficients[ ,"Value"]))

# combine with coefficient and p_value
printCoefmat(coefficients[ ,c("Value", "Std. Error", "odds_ratio", "p_value")], 
             P.values=TRUE, 
             has.Pvalue=TRUE)
##                           Value  Std. Error odds_ratio   p_value    
## A.Mismatch.No.       0.00050392  0.01055601     1.0005  0.961925    
## B.Mismatch.No.       0.01370436  0.01208611     1.0138  0.256839    
## C.Mismatch.No.      -0.02600325  0.01449154     0.9743  0.072753 .  
## DRB1.Mismatch.No.    0.03668191  0.01308992     1.0374  0.005074 ** 
## DRB345.Mismatch.No. -0.03438792  0.00964633     0.9662  0.000364 ***
## DP.Mismatch.No.     -0.02142974  0.00948431     0.9788  0.023853 *  
## DQ.Mismatch.No.      0.03782805  0.01135264     1.0386  0.000862 ***
## age                 -0.01709165  0.00609452     0.9831  0.005041 ** 
## Time_biopsy         -0.00219586  0.00458371     0.9978  0.631898    
## 1|2                 -0.95908914  0.42037203     0.3832  0.022517 *  
## 2|3                  2.34344035  0.43440975    10.4170 6.870e-08 ***
## 3|4                  6.30458471  1.08049763   547.0743 5.383e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
library(brant)
brant(model)
## Warning: glm.fit: algorithm did not converge
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## ---------------------------------------------------- 
## Test for     X2  df  probability 
## ---------------------------------------------------- 
## Omnibus          8.3 18  0.97
## A.Mismatch.No.       0   2   1
## B.Mismatch.No.       3.06    2   0.22
## C.Mismatch.No.       0   2   1
## DRB1.Mismatch.No.    -0.12   2   1
## DRB345.Mismatch.No.  0.7 2   0.71
## DP.Mismatch.No.  -0.01   2   1
## DQ.Mismatch.No.  0.06    2   0.97
## age          0   2   1
## Time_biopsy      0   2   1
## ---------------------------------------------------- 
## 
## H0: Parallel Regression Assumption holds
df_long$trans_id<-df_long%>%filter()
glm(pom_formula, data = df_long, family = binomial(link = "logit"))
## 
## Call:  glm(formula = pom_formula, family = binomial(link = "logit"), 
##     data = df_long)
## 
## Coefficients:
##         (Intercept)       A.Mismatch.No.       B.Mismatch.No.  
##            0.895554            -0.001885             0.009735  
##      C.Mismatch.No.    DRB1.Mismatch.No.  DRB345.Mismatch.No.  
##           -0.028640             0.036148            -0.036378  
##     DP.Mismatch.No.      DQ.Mismatch.No.                  age  
##           -0.016782             0.039024            -0.015641  
##         Time_biopsy  
##            0.001610  
## 
## Degrees of Freedom: 897 Total (i.e. Null);  888 Residual
## Null Deviance:       1219 
## Residual Deviance: 1178  AIC: 1198
#repolr(pom_formula, subjects = "trans_id",data = df_long, categories = 4, times = c(1, 2, 3, 4, 6, 8, 10, 12, 26, 52), corstr = "ar1", 
#        tol = 0.001, scalevalue = 1, alpha = 0.5, po.test=TRUE, fixed=FALSE) 

po with max grades

df$Max_Grade<-df$Max.Grade %>% as.factor()
pom_formula2 <- create_pom_formula("Max_Grade", independent_vars, "time2fail" ,df)
model2 <- polr(pom_formula2, data = df, Hess = TRUE)
coefficients <- summary(model2)$coefficients
p_value <- (1 - pnorm(abs(coefficients[ ,"t value"]), 0, 1))*2
coefficients <- cbind(coefficients, p_value)
coefficients <- cbind(coefficients, odds_ratio = exp(coefficients[ ,"Value"]))
printCoefmat(coefficients[ ,c("Value", "Std. Error", "odds_ratio", "p_value")], 
             P.values=TRUE, 
             has.Pvalue=TRUE)
##                          Value Std. Error odds_ratio   p_value    
## A.Mismatch.No.       0.0180628  0.0362504     1.0182 0.6182876    
## B.Mismatch.No.       0.0313153  0.0399608     1.0318 0.4332454    
## C.Mismatch.No.      -0.0084999  0.0498731     0.9915 0.8646724    
## DRB1.Mismatch.No.    0.0444439  0.0456208     1.0454 0.3299557    
## DRB345.Mismatch.No. -0.0282089  0.0319086     0.9722 0.3766672    
## DP.Mismatch.No.     -0.0713978  0.0334554     0.9311 0.0328330 *  
## DQ.Mismatch.No.      0.0840075  0.0383948     1.0876 0.0286700 *  
## age                 -0.0567930  0.0216237     0.9448 0.0086286 ** 
## 0|1                 -6.2233976  1.6339131     0.0020 0.0001396 ***
## 1|2                 -1.0775802  1.4072995     0.3404 0.4438502    
## 2|3                  3.3304985  1.6404126    27.9523 0.0423280 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
conf<-exp(confint(model2, level = 0.90))
## Waiting for profiling to be done...
var<-c("A.Mismatch.No.", "B.Mismatch.No.", "C.Mismatch.No.", "", "DP.Mismatch.No.", "DQ.Mismatch.No.", "age")
lower<-c(conf[,1])
upper<-c(conf[,2])
odds<-exp(coefficients[ 0:8,"Value"])
ndf <- data.frame(Variable = independent_vars, Lower = lower, Upper = upper, Odds = odds)
ndf$err <- ndf$Upper - ndf$Odds

svg("figs/forest_plot.svg")
ndf %>% forestplot(labeltext=Variable, mean=Odds, lower=Lower, upper=Upper, xlab="Odds Ratio", zero = 1,
             cex  = 2, lineheight = "auto")
dev.off()
## png 
##   2
df$X1st.DSA<-anydate(df$X1st.DSA)
df$X1st.DSA_status<-ifelse(df$X1st.DSA==anydate("2024-03-10"), 0, 1)
df$X1st.DSA<-ifelse(is.na(df$X1st.DSA), anydate("2024-03-10"), df$X1st.DSA)
df$X1st.DSA_time<-anydate(df$X1st.DSA)-df$Transplant.Date
df$X1st.DSA_status
##  [1]  1 NA NA NA NA NA NA NA NA NA NA NA NA NA NA  1  1 NA  1 NA NA NA NA NA NA
## [26] NA  1 NA NA NA NA NA NA NA NA NA NA NA NA NA  1 NA  1 NA  1 NA  1 NA NA NA
## [51] NA NA NA NA  1 NA NA NA NA  1 NA NA  1 NA NA  1 NA NA NA NA NA NA NA NA NA
## [76]  1 NA NA NA NA  1 NA  1  1 NA NA NA NA NA NA NA NA NA NA NA NA NA NA
df$Transplant.Date
##  [1] "2023-06-23" "2023-03-22" "2023-02-20" "2022-10-16" "2022-07-24"
##  [6] "2022-04-10" "2022-03-20" "2022-02-02" "2022-01-09" "2022-01-07"
## [11] "2021-11-07" "2021-09-09" "2021-09-02" "2021-07-16" "2021-06-07"
## [16] "2021-04-19" "2020-11-16" "2020-11-07" "2020-09-15" "2020-08-20"
## [21] "2020-07-04" "2020-06-27" "2020-04-13" "2020-03-11" "2020-01-16"
## [26] "2020-01-02" "2019-11-04" "2019-10-24" "2019-09-30" "2019-09-28"
## [31] "2019-09-11" "2019-09-01" "2019-08-14" "2019-07-29" "2019-06-20"
## [36] "2019-06-16" "2019-04-07" "2019-02-16" "2019-01-05" "2018-12-19"
## [41] "2018-12-11" "2018-12-04" "2018-07-19" "2018-06-29" "2018-06-28"
## [46] "2018-06-28" "2018-06-22" "2018-04-04" "2018-02-09" "2018-02-04"
## [51] "2017-11-08" "2017-10-29" "2017-10-29" "2017-10-09" "2017-09-15"
## [56] "2017-08-15" "2017-05-24" "2017-05-09" "2017-03-14" "2017-02-26"
## [61] "2017-02-22" "2017-02-18" "2016-11-06" "2016-10-17" "2016-09-12"
## [66] "2016-07-31" "2016-05-15" "2016-05-07" "2016-04-26" "2016-03-27"
## [71] "2015-12-11" "2015-09-25" "2015-09-14" "2015-06-26" "2015-04-02"
## [76] "2015-03-12" "2015-02-02" "2015-01-31" "2014-09-30" "2014-09-15"
## [81] "2014-08-15" "2014-08-02" "2014-06-19" "2014-03-26" "2013-11-29"
## [86] "2013-10-22" "2013-08-08" "2013-07-07" "2013-02-23" "2012-07-14"
## [91] "2012-06-30" "2011-10-21" "2011-10-18" "2011-08-19" "2011-07-01"
## [96] "2011-06-23" "2011-02-10" "2010-11-15"
ggsurvplot(surv_fit(Surv(X1st.DSA_time, X1st.DSA_status) ~ 1, data = df), 
           data = df, 
           xlab = "Time to Failure", 
           ylab = "Survival Probability",
           title = "Kaplan-Meier Survival Curves",)

univariate_coxph(df, independent_vars, "X1st.DSA_status", "X1st.DSA_time")

## 
## Cox proportional hazards model for A.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
##                   coef exp(coef) se(coef)     z Pr(>|z|)
## A.Mismatch.No. 0.01999   1.02019  0.03804 0.525    0.599
## 
##                exp(coef) exp(-coef) lower .95 upper .95
## A.Mismatch.No.      1.02     0.9802    0.9469     1.099
## 
## Concordance= 0.5  (se = 0.077 )
## Likelihood ratio test= 0.27  on 1 df,   p=0.6
## Wald test            = 0.28  on 1 df,   p=0.6
## Score (logrank) test = 0.28  on 1 df,   p=0.6

## 
## Cox proportional hazards model for B.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
##                   coef exp(coef) se(coef)     z Pr(>|z|)
## B.Mismatch.No. 0.04676   1.04787  0.04196 1.114    0.265
## 
##                exp(coef) exp(-coef) lower .95 upper .95
## B.Mismatch.No.     1.048     0.9543    0.9651     1.138
## 
## Concordance= 0.665  (se = 0.071 )
## Likelihood ratio test= 1.26  on 1 df,   p=0.3
## Wald test            = 1.24  on 1 df,   p=0.3
## Score (logrank) test = 1.25  on 1 df,   p=0.3

## 
## Cox proportional hazards model for C.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
##                    coef exp(coef) se(coef)      z Pr(>|z|)
## C.Mismatch.No. -0.03032   0.97014  0.05255 -0.577    0.564
## 
##                exp(coef) exp(-coef) lower .95 upper .95
## C.Mismatch.No.    0.9701      1.031    0.8752     1.075
## 
## Concordance= 0.515  (se = 0.071 )
## Likelihood ratio test= 0.33  on 1 df,   p=0.6
## Wald test            = 0.33  on 1 df,   p=0.6
## Score (logrank) test = 0.33  on 1 df,   p=0.6

## 
## Cox proportional hazards model for DRB1.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
##                      coef exp(coef) se(coef)     z Pr(>|z|)
## DRB1.Mismatch.No. 0.02296   1.02323  0.04633 0.496     0.62
## 
##                   exp(coef) exp(-coef) lower .95 upper .95
## DRB1.Mismatch.No.     1.023     0.9773    0.9344      1.12
## 
## Concordance= 0.548  (se = 0.111 )
## Likelihood ratio test= 0.25  on 1 df,   p=0.6
## Wald test            = 0.25  on 1 df,   p=0.6
## Score (logrank) test = 0.25  on 1 df,   p=0.6

## 
## Cox proportional hazards model for DRB345.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
##                        coef exp(coef) se(coef)     z Pr(>|z|)
## DRB345.Mismatch.No. 0.01726   1.01741  0.03941 0.438    0.661
## 
##                     exp(coef) exp(-coef) lower .95 upper .95
## DRB345.Mismatch.No.     1.017     0.9829    0.9418     1.099
## 
## Concordance= 0.529  (se = 0.083 )
## Likelihood ratio test= 0.19  on 1 df,   p=0.7
## Wald test            = 0.19  on 1 df,   p=0.7
## Score (logrank) test = 0.19  on 1 df,   p=0.7

## 
## Cox proportional hazards model for DP.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
##                     coef exp(coef) se(coef)      z Pr(>|z|)
## DP.Mismatch.No. -0.01717   0.98298  0.03961 -0.434    0.665
## 
##                 exp(coef) exp(-coef) lower .95 upper .95
## DP.Mismatch.No.     0.983      1.017    0.9096     1.062
## 
## Concordance= 0.485  (se = 0.082 )
## Likelihood ratio test= 0.19  on 1 df,   p=0.7
## Wald test            = 0.19  on 1 df,   p=0.7
## Score (logrank) test = 0.19  on 1 df,   p=0.7

## 
## Cox proportional hazards model for DQ.Mismatch.No. :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
##                     coef exp(coef) se(coef)      z Pr(>|z|)
## DQ.Mismatch.No. -0.02360   0.97668  0.03926 -0.601    0.548
## 
##                 exp(coef) exp(-coef) lower .95 upper .95
## DQ.Mismatch.No.    0.9767      1.024    0.9044     1.055
## 
## Concordance= 0.562  (se = 0.097 )
## Likelihood ratio test= 0.36  on 1 df,   p=0.5
## Wald test            = 0.36  on 1 df,   p=0.5
## Score (logrank) test = 0.36  on 1 df,   p=0.5

## 
## Cox proportional hazards model for age :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
##         coef exp(coef) se(coef)     z Pr(>|z|)
## age -0.01713   0.98302  0.02255 -0.76    0.448
## 
##     exp(coef) exp(-coef) lower .95 upper .95
## age     0.983      1.017    0.9405     1.027
## 
## Concordance= 0.632  (se = 0.088 )
## Likelihood ratio test= 0.52  on 1 df,   p=0.5
## Wald test            = 0.58  on 1 df,   p=0.4
## Score (logrank) test = 0.58  on 1 df,   p=0.4
## 
##                                Variable        HR  Lower_CI Upper_CI   p_value
## A.Mismatch.No.           A.Mismatch.No. 1.0201888 0.9468937 1.099157 0.5992739
## B.Mismatch.No.           B.Mismatch.No. 1.0478670 0.9651414 1.137683 0.2651276
## C.Mismatch.No.           C.Mismatch.No. 0.9701363 0.8751862 1.075388 0.5639887
## DRB1.Mismatch.No.     DRB1.Mismatch.No. 1.0232300 0.9344066 1.120497 0.6201388
## DRB345.Mismatch.No. DRB345.Mismatch.No. 1.0174110 0.9417806 1.099115 0.6614020
## DP.Mismatch.No.         DP.Mismatch.No. 0.9829767 0.9095570 1.062323 0.6646446
## DQ.Mismatch.No.         DQ.Mismatch.No. 0.9766793 0.9043500 1.054793 0.5477780
## age                                 age 0.9830200 0.9405224 1.027438 0.4475451
##                       C_index
## A.Mismatch.No.      0.5000000
## B.Mismatch.No.      0.6654412
## C.Mismatch.No.      0.5147059
## DRB1.Mismatch.No.   0.5477941
## DRB345.Mismatch.No. 0.5294118
## DP.Mismatch.No.     0.4852941
## DQ.Mismatch.No.     0.5625000
## age                 0.6323529
## $A.Mismatch.No.
## $A.Mismatch.No.$Model
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##                   coef exp(coef) se(coef)     z     p
## A.Mismatch.No. 0.01999   1.02019  0.03804 0.525 0.599
## 
## Likelihood ratio test=0.27  on 1 df, p=0.6012
## n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
## $A.Mismatch.No.$Summary
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
##                   coef exp(coef) se(coef)     z Pr(>|z|)
## A.Mismatch.No. 0.01999   1.02019  0.03804 0.525    0.599
## 
##                exp(coef) exp(-coef) lower .95 upper .95
## A.Mismatch.No.      1.02     0.9802    0.9469     1.099
## 
## Concordance= 0.5  (se = 0.077 )
## Likelihood ratio test= 0.27  on 1 df,   p=0.6
## Wald test            = 0.28  on 1 df,   p=0.6
## Score (logrank) test = 0.28  on 1 df,   p=0.6
## 
## 
## 
## $B.Mismatch.No.
## $B.Mismatch.No.$Model
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##                   coef exp(coef) se(coef)     z     p
## B.Mismatch.No. 0.04676   1.04787  0.04196 1.114 0.265
## 
## Likelihood ratio test=1.26  on 1 df, p=0.2607
## n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
## $B.Mismatch.No.$Summary
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
##                   coef exp(coef) se(coef)     z Pr(>|z|)
## B.Mismatch.No. 0.04676   1.04787  0.04196 1.114    0.265
## 
##                exp(coef) exp(-coef) lower .95 upper .95
## B.Mismatch.No.     1.048     0.9543    0.9651     1.138
## 
## Concordance= 0.665  (se = 0.071 )
## Likelihood ratio test= 1.26  on 1 df,   p=0.3
## Wald test            = 1.24  on 1 df,   p=0.3
## Score (logrank) test = 1.25  on 1 df,   p=0.3
## 
## 
## 
## $C.Mismatch.No.
## $C.Mismatch.No.$Model
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##                    coef exp(coef) se(coef)      z     p
## C.Mismatch.No. -0.03032   0.97014  0.05255 -0.577 0.564
## 
## Likelihood ratio test=0.33  on 1 df, p=0.5642
## n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
## $C.Mismatch.No.$Summary
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
##                    coef exp(coef) se(coef)      z Pr(>|z|)
## C.Mismatch.No. -0.03032   0.97014  0.05255 -0.577    0.564
## 
##                exp(coef) exp(-coef) lower .95 upper .95
## C.Mismatch.No.    0.9701      1.031    0.8752     1.075
## 
## Concordance= 0.515  (se = 0.071 )
## Likelihood ratio test= 0.33  on 1 df,   p=0.6
## Wald test            = 0.33  on 1 df,   p=0.6
## Score (logrank) test = 0.33  on 1 df,   p=0.6
## 
## 
## 
## $DRB1.Mismatch.No.
## $DRB1.Mismatch.No.$Model
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##                      coef exp(coef) se(coef)     z    p
## DRB1.Mismatch.No. 0.02296   1.02323  0.04633 0.496 0.62
## 
## Likelihood ratio test=0.25  on 1 df, p=0.6156
## n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
## $DRB1.Mismatch.No.$Summary
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
##                      coef exp(coef) se(coef)     z Pr(>|z|)
## DRB1.Mismatch.No. 0.02296   1.02323  0.04633 0.496     0.62
## 
##                   exp(coef) exp(-coef) lower .95 upper .95
## DRB1.Mismatch.No.     1.023     0.9773    0.9344      1.12
## 
## Concordance= 0.548  (se = 0.111 )
## Likelihood ratio test= 0.25  on 1 df,   p=0.6
## Wald test            = 0.25  on 1 df,   p=0.6
## Score (logrank) test = 0.25  on 1 df,   p=0.6
## 
## 
## 
## $DRB345.Mismatch.No.
## $DRB345.Mismatch.No.$Model
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##                        coef exp(coef) se(coef)     z     p
## DRB345.Mismatch.No. 0.01726   1.01741  0.03941 0.438 0.661
## 
## Likelihood ratio test=0.19  on 1 df, p=0.6607
## n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
## $DRB345.Mismatch.No.$Summary
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
##                        coef exp(coef) se(coef)     z Pr(>|z|)
## DRB345.Mismatch.No. 0.01726   1.01741  0.03941 0.438    0.661
## 
##                     exp(coef) exp(-coef) lower .95 upper .95
## DRB345.Mismatch.No.     1.017     0.9829    0.9418     1.099
## 
## Concordance= 0.529  (se = 0.083 )
## Likelihood ratio test= 0.19  on 1 df,   p=0.7
## Wald test            = 0.19  on 1 df,   p=0.7
## Score (logrank) test = 0.19  on 1 df,   p=0.7
## 
## 
## 
## $DP.Mismatch.No.
## $DP.Mismatch.No.$Model
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##                     coef exp(coef) se(coef)      z     p
## DP.Mismatch.No. -0.01717   0.98298  0.03961 -0.434 0.665
## 
## Likelihood ratio test=0.19  on 1 df, p=0.6648
## n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
## $DP.Mismatch.No.$Summary
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
##                     coef exp(coef) se(coef)      z Pr(>|z|)
## DP.Mismatch.No. -0.01717   0.98298  0.03961 -0.434    0.665
## 
##                 exp(coef) exp(-coef) lower .95 upper .95
## DP.Mismatch.No.     0.983      1.017    0.9096     1.062
## 
## Concordance= 0.485  (se = 0.082 )
## Likelihood ratio test= 0.19  on 1 df,   p=0.7
## Wald test            = 0.19  on 1 df,   p=0.7
## Score (logrank) test = 0.19  on 1 df,   p=0.7
## 
## 
## 
## $DQ.Mismatch.No.
## $DQ.Mismatch.No.$Model
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##                     coef exp(coef) se(coef)      z     p
## DQ.Mismatch.No. -0.02360   0.97668  0.03926 -0.601 0.548
## 
## Likelihood ratio test=0.36  on 1 df, p=0.5479
## n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
## $DQ.Mismatch.No.$Summary
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
##                     coef exp(coef) se(coef)      z Pr(>|z|)
## DQ.Mismatch.No. -0.02360   0.97668  0.03926 -0.601    0.548
## 
##                 exp(coef) exp(-coef) lower .95 upper .95
## DQ.Mismatch.No.    0.9767      1.024    0.9044     1.055
## 
## Concordance= 0.562  (se = 0.097 )
## Likelihood ratio test= 0.36  on 1 df,   p=0.5
## Wald test            = 0.36  on 1 df,   p=0.5
## Score (logrank) test = 0.36  on 1 df,   p=0.5
## 
## 
## 
## $age
## $age$Model
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##         coef exp(coef) se(coef)     z     p
## age -0.01713   0.98302  0.02255 -0.76 0.448
## 
## Likelihood ratio test=0.52  on 1 df, p=0.4694
## n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
## $age$Summary
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
##         coef exp(coef) se(coef)     z Pr(>|z|)
## age -0.01713   0.98302  0.02255 -0.76    0.448
## 
##     exp(coef) exp(-coef) lower .95 upper .95
## age     0.983      1.017    0.9405     1.027
## 
## Concordance= 0.632  (se = 0.088 )
## Likelihood ratio test= 0.52  on 1 df,   p=0.5
## Wald test            = 0.58  on 1 df,   p=0.4
## Score (logrank) test = 0.58  on 1 df,   p=0.4
plot_quartiles_and_logrank(df, independent_vars, "X1st.DSA_status", "X1st.DSA_time")
## Log-rank test result for A.Mismatch.No. :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
## n=17, 81 observations deleted due to missingness.
## 
##                  N Observed Expected (O-E)^2/E (O-E)^2/V
## A.Mismatch.No.=1 9        9     9.39    0.0165    0.0422
## A.Mismatch.No.=2 8        8     7.61    0.0204    0.0422
## 
##  Chisq= 0  on 1 degrees of freedom, p= 0.8

## Log-rank test result for B.Mismatch.No. :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
## n=17, 81 observations deleted due to missingness.
## 
##                  N Observed Expected (O-E)^2/E (O-E)^2/V
## B.Mismatch.No.=1 8        8     9.93     0.374     0.955
## B.Mismatch.No.=2 9        9     7.07     0.526     0.955
## 
##  Chisq= 1  on 1 degrees of freedom, p= 0.3

## Log-rank test result for C.Mismatch.No. :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
## n=17, 81 observations deleted due to missingness.
## 
##                  N Observed Expected (O-E)^2/E (O-E)^2/V
## C.Mismatch.No.=1 8        8      5.1     1.652      2.79
## C.Mismatch.No.=2 9        9     11.9     0.708      2.79
## 
##  Chisq= 2.8  on 1 degrees of freedom, p= 0.1

## Log-rank test result for DRB1.Mismatch.No. :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
## n=17, 81 observations deleted due to missingness.
## 
##                     N Observed Expected (O-E)^2/E (O-E)^2/V
## DRB1.Mismatch.No.=1 9        9     8.81   0.00405   0.00917
## DRB1.Mismatch.No.=2 8        8     8.19   0.00436   0.00917
## 
##  Chisq= 0  on 1 degrees of freedom, p= 0.9

## Log-rank test result for DRB345.Mismatch.No. :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
## n=17, 81 observations deleted due to missingness.
## 
##                       N Observed Expected (O-E)^2/E (O-E)^2/V
## DRB345.Mismatch.No.=1 9        9     8.53    0.0263    0.0582
## DRB345.Mismatch.No.=2 8        8     8.47    0.0265    0.0582
## 
##  Chisq= 0.1  on 1 degrees of freedom, p= 0.8

## Log-rank test result for DP.Mismatch.No. :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
## n=17, 81 observations deleted due to missingness.
## 
##                   N Observed Expected (O-E)^2/E (O-E)^2/V
## DP.Mismatch.No.=1 8        8     6.61     0.292     0.536
## DP.Mismatch.No.=2 9        9    10.39     0.186     0.536
## 
##  Chisq= 0.5  on 1 degrees of freedom, p= 0.5

## Log-rank test result for DQ.Mismatch.No. :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
## n=17, 81 observations deleted due to missingness.
## 
##                    N Observed Expected (O-E)^2/E (O-E)^2/V
## DQ.Mismatch.No.=1  6        6     5.71   0.01479    0.0235
## DQ.Mismatch.No.=2 11       11    11.29   0.00748    0.0235
## 
##  Chisq= 0  on 1 degrees of freedom, p= 0.9

## Log-rank test result for age :
## Call:
## survdiff(formula = fit_formula, data = data)
## 
## n=17, 81 observations deleted due to missingness.
## 
##        N Observed Expected (O-E)^2/E (O-E)^2/V
## age=1 10       10     8.23     0.380     0.807
## age=2  7        7     8.77     0.357     0.807
## 
##  Chisq= 0.8  on 1 degrees of freedom, p= 0.4

## $A.Mismatch.No.
## $A.Mismatch.No.$plot

## 
## $A.Mismatch.No.$log_rank_result
## Call:
## survdiff(formula = fit_formula, data = data)
## 
## n=17, 81 observations deleted due to missingness.
## 
##                  N Observed Expected (O-E)^2/E (O-E)^2/V
## A.Mismatch.No.=1 9        9     9.39    0.0165    0.0422
## A.Mismatch.No.=2 8        8     7.61    0.0204    0.0422
## 
##  Chisq= 0  on 1 degrees of freedom, p= 0.8 
## 
## 
## $B.Mismatch.No.
## $B.Mismatch.No.$plot

## 
## $B.Mismatch.No.$log_rank_result
## Call:
## survdiff(formula = fit_formula, data = data)
## 
## n=17, 81 observations deleted due to missingness.
## 
##                  N Observed Expected (O-E)^2/E (O-E)^2/V
## B.Mismatch.No.=1 8        8     9.93     0.374     0.955
## B.Mismatch.No.=2 9        9     7.07     0.526     0.955
## 
##  Chisq= 1  on 1 degrees of freedom, p= 0.3 
## 
## 
## $C.Mismatch.No.
## $C.Mismatch.No.$plot

## 
## $C.Mismatch.No.$log_rank_result
## Call:
## survdiff(formula = fit_formula, data = data)
## 
## n=17, 81 observations deleted due to missingness.
## 
##                  N Observed Expected (O-E)^2/E (O-E)^2/V
## C.Mismatch.No.=1 8        8      5.1     1.652      2.79
## C.Mismatch.No.=2 9        9     11.9     0.708      2.79
## 
##  Chisq= 2.8  on 1 degrees of freedom, p= 0.1 
## 
## 
## $DRB1.Mismatch.No.
## $DRB1.Mismatch.No.$plot

## 
## $DRB1.Mismatch.No.$log_rank_result
## Call:
## survdiff(formula = fit_formula, data = data)
## 
## n=17, 81 observations deleted due to missingness.
## 
##                     N Observed Expected (O-E)^2/E (O-E)^2/V
## DRB1.Mismatch.No.=1 9        9     8.81   0.00405   0.00917
## DRB1.Mismatch.No.=2 8        8     8.19   0.00436   0.00917
## 
##  Chisq= 0  on 1 degrees of freedom, p= 0.9 
## 
## 
## $DRB345.Mismatch.No.
## $DRB345.Mismatch.No.$plot

## 
## $DRB345.Mismatch.No.$log_rank_result
## Call:
## survdiff(formula = fit_formula, data = data)
## 
## n=17, 81 observations deleted due to missingness.
## 
##                       N Observed Expected (O-E)^2/E (O-E)^2/V
## DRB345.Mismatch.No.=1 9        9     8.53    0.0263    0.0582
## DRB345.Mismatch.No.=2 8        8     8.47    0.0265    0.0582
## 
##  Chisq= 0.1  on 1 degrees of freedom, p= 0.8 
## 
## 
## $DP.Mismatch.No.
## $DP.Mismatch.No.$plot

## 
## $DP.Mismatch.No.$log_rank_result
## Call:
## survdiff(formula = fit_formula, data = data)
## 
## n=17, 81 observations deleted due to missingness.
## 
##                   N Observed Expected (O-E)^2/E (O-E)^2/V
## DP.Mismatch.No.=1 8        8     6.61     0.292     0.536
## DP.Mismatch.No.=2 9        9    10.39     0.186     0.536
## 
##  Chisq= 0.5  on 1 degrees of freedom, p= 0.5 
## 
## 
## $DQ.Mismatch.No.
## $DQ.Mismatch.No.$plot

## 
## $DQ.Mismatch.No.$log_rank_result
## Call:
## survdiff(formula = fit_formula, data = data)
## 
## n=17, 81 observations deleted due to missingness.
## 
##                    N Observed Expected (O-E)^2/E (O-E)^2/V
## DQ.Mismatch.No.=1  6        6     5.71   0.01479    0.0235
## DQ.Mismatch.No.=2 11       11    11.29   0.00748    0.0235
## 
##  Chisq= 0  on 1 degrees of freedom, p= 0.9 
## 
## 
## $age
## $age$plot

## 
## $age$log_rank_result
## Call:
## survdiff(formula = fit_formula, data = data)
## 
## n=17, 81 observations deleted due to missingness.
## 
##        N Observed Expected (O-E)^2/E (O-E)^2/V
## age=1 10       10     8.23     0.380     0.807
## age=2  7        7     8.77     0.357     0.807
## 
##  Chisq= 0.8  on 1 degrees of freedom, p= 0.4

plot density curves

plot_density_faceted <- function(data, colnames) {
  # Check if all colnames exist in the dataframe
  missing_cols <- colnames[!colnames %in% names(data)]
  if(length(missing_cols) > 0) {
    stop("The following columns are missing from the dataframe: ", paste(missing_cols, collapse = ", "), ".")
  }
  
  # Melt the dataframe to long format for faceting
  long_data <- reshape2::melt(data, measure.vars = colnames)
  
  # Create the density plot
  p <- ggplot(long_data, aes(x = value)) +
    geom_density(aes(fill = variable), alpha = 0.5) +
    facet_wrap(~ variable, scales = "free") +
    labs(x = "Value", y = "Density") +
    theme_minimal()
  
  svg("figs/density_plot.svg")
  print(p)
  dev.off()
}
plot_density_faceted(df, HLA_vars)
## png 
##   2
plot_density_faceted(df, HLA_al_vars)
## png 
##   2

DSA survival curve (eplet)

cox_results<-univariate_coxph(df, independent_al_vars, "X1st.DSA_status", "X1st.DSA_time")

## 
## Cox proportional hazards model for A :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
##     coef exp(coef) se(coef)     z Pr(>|z|)
## A 0.9612    2.6149   0.5916 1.625    0.104
## 
##   exp(coef) exp(-coef) lower .95 upper .95
## A     2.615     0.3824    0.8201     8.337
## 
## Concordance= 0.618  (se = 0.063 )
## Likelihood ratio test= 2.96  on 1 df,   p=0.09
## Wald test            = 2.64  on 1 df,   p=0.1
## Score (logrank) test = 2.83  on 1 df,   p=0.09

## 
## Cox proportional hazards model for B :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
##     coef exp(coef) se(coef)    z Pr(>|z|)
## B 0.1525    1.1648   0.6644 0.23    0.818
## 
##   exp(coef) exp(-coef) lower .95 upper .95
## B     1.165     0.8585    0.3168     4.283
## 
## Concordance= 0.559  (se = 0.044 )
## Likelihood ratio test= 0.05  on 1 df,   p=0.8
## Wald test            = 0.05  on 1 df,   p=0.8
## Score (logrank) test = 0.05  on 1 df,   p=0.8

## 
## Cox proportional hazards model for C :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
##      coef exp(coef) se(coef)      z Pr(>|z|)
## C -0.3245    0.7229   0.3482 -0.932    0.351
## 
##   exp(coef) exp(-coef) lower .95 upper .95
## C    0.7229      1.383    0.3653      1.43
## 
## Concordance= 0.562  (se = 0.084 )
## Likelihood ratio test= 0.83  on 1 df,   p=0.4
## Wald test            = 0.87  on 1 df,   p=0.4
## Score (logrank) test = 0.88  on 1 df,   p=0.3

## 
## Cox proportional hazards model for DRB1 :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
##          coef exp(coef) se(coef)      z Pr(>|z|)
## DRB1 -0.06552   0.93658  0.42950 -0.153    0.879
## 
##      exp(coef) exp(-coef) lower .95 upper .95
## DRB1    0.9366      1.068    0.4036     2.173
## 
## Concordance= 0.485  (se = 0.062 )
## Likelihood ratio test= 0.02  on 1 df,   p=0.9
## Wald test            = 0.02  on 1 df,   p=0.9
## Score (logrank) test = 0.02  on 1 df,   p=0.9

## 
## Cox proportional hazards model for DRB345 :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
##           coef exp(coef) se(coef)      z Pr(>|z|)
## DRB345 -0.1247    0.8827   0.1385 -0.901    0.368
## 
##        exp(coef) exp(-coef) lower .95 upper .95
## DRB345    0.8827      1.133    0.6729     1.158
## 
## Concordance= 0.57  (se = 0.08 )
## Likelihood ratio test= 0.82  on 1 df,   p=0.4
## Wald test            = 0.81  on 1 df,   p=0.4
## Score (logrank) test = 0.82  on 1 df,   p=0.4

## 
## Cox proportional hazards model for DP :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
##       coef exp(coef) se(coef)      z Pr(>|z|)  
## DP -0.6149    0.5407   0.3030 -2.029   0.0424 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##    exp(coef) exp(-coef) lower .95 upper .95
## DP    0.5407      1.849    0.2986    0.9793
## 
## Concordance= 0.614  (se = 0.094 )
## Likelihood ratio test= 4.64  on 1 df,   p=0.03
## Wald test            = 4.12  on 1 df,   p=0.04
## Score (logrank) test = 4.38  on 1 df,   p=0.04

## 
## Cox proportional hazards model for DQ :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
##       coef exp(coef) se(coef)      z Pr(>|z|)
## DQ -0.2390    0.7874   0.2707 -0.883    0.377
## 
##    exp(coef) exp(-coef) lower .95 upper .95
## DQ    0.7874       1.27    0.4633     1.338
## 
## Concordance= 0.588  (se = 0.082 )
## Likelihood ratio test= 0.76  on 1 df,   p=0.4
## Wald test            = 0.78  on 1 df,   p=0.4
## Score (logrank) test = 0.79  on 1 df,   p=0.4

## 
## Cox proportional hazards model for age :
## Call:
## coxph(formula = fit_formula, data = data, na.action = na.exclude)
## 
##   n= 17, number of events= 17 
##    (81 observations deleted due to missingness)
## 
##         coef exp(coef) se(coef)     z Pr(>|z|)
## age -0.01713   0.98302  0.02255 -0.76    0.448
## 
##     exp(coef) exp(-coef) lower .95 upper .95
## age     0.983      1.017    0.9405     1.027
## 
## Concordance= 0.632  (se = 0.088 )
## Likelihood ratio test= 0.52  on 1 df,   p=0.5
## Wald test            = 0.58  on 1 df,   p=0.4
## Score (logrank) test = 0.58  on 1 df,   p=0.4
## 
##        Variable        HR  Lower_CI  Upper_CI    p_value   C_index
## A             A 2.6148671 0.8201377 8.3370509 0.10420631 0.6176471
## B             B 1.1647798 0.3167629 4.2830516 0.81841013 0.5588235
## C             C 0.7228977 0.3653291 1.4304392 0.35139506 0.5625000
## DRB1       DRB1 0.9365844 0.4036056 2.1733848 0.87876111 0.4852941
## DRB345   DRB345 0.8827405 0.6729353 1.1579580 0.36770740 0.5698529
## DP           DP 0.5407217 0.2985717 0.9792619 0.04244575 0.6139706
## DQ           DQ 0.7874395 0.4632733 1.3384345 0.37727129 0.5882353
## age         age 0.9830200 0.9405224 1.0274379 0.44754506 0.6323529